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CENTER  FOR  ANALYSIS  OF  HETEROGENEOUS  AND  NONLINEAR  MEDIA 


I.  Introduction 

/The  Center  for  Analysis  of  Heterogeneous  and  Nonlinear  Media  was  established  in  October  1986 
through  the  University  Research  Initiative  (URI)  Program.  Funding  for  the  Center  is  primarily  from  the 
Air  Force  Office  of  Scientific  Research  through  the  URI  Program  and  the  Applied  Analysis  Program 
managed  by  Atje  Nachman. 

The  Center  consists  of  a  number  of  faculty,  postdoctoral  researchers,  visiting  scientists  and  graduate 
students  integrated  within  the  Courant  Institute.  Establishment  of  this  Center  provided  a  scientific  focus 
and  a  unified  computational  environment  of  Sun  workstations.  The  result  has  been  a  very  active  and  col¬ 
laborative  research  group.  Numerical  computation  has  played  an  extremely  important  role  in  the  research 
of  this  Center,  for  simulation  of  physical  systems  and  as  an  inspiration  and  guide  for  modeling  and 
analysis. 

The  Courant  Institute  has  provided  a  fruitful  setting  for  the  Center.  Many  faculty,  postdoctoral  visi¬ 
tors  and  graduate  students  of  the  Institute  have  become  involved  in  the  research  projects  of  the  Center. 
The  Institute’s  extensive  system  of  seminars  and  courses  has  helped  keep  the  Center  at  the  forefront  of 
current  scientific  developments. 

The  scientific  research  of  the  Center  has  consisted  of  modeling,  mathematical  analysis  and  compu¬ 
tation  of  problems  in  composite  and  random  materials,  nonlinear  optics  and  fluid  dynamics.  The  results 
of  this  research  are  presented  in  a  number  of  research  publications  and  several  books.  Two  workshops,  on 
vortex  dynamics  and  on  composite  materials,  were  organized  by  the  Center.  These  workshops  provided  a 
very  useful  forum  for  presentation  of  the  results  of  the  Center,  as  well  as  for  scientific  interchange  among 
a  wide  range  of  participants. 
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II.  Focus  and  Goals 

\  The  focus  of  the  Center  has  been  on  the  following  research  areas: 

•  Vortex  dynamics  in  incompressible  fluids 

•  Numerical  methods  for  fluid  dynamics 

•  Macroscopic  properties  of  composite  materials 

•  Wave  propagation  in  random  media 

•  Nonlinear  phenomena  in  optical  materials.  '  f  , 

Basic  research  in  these  subjects  is  of  great  technological  importance  and  poses  mathematical,  numerical 
and  physical  challenges.  In  the  research  of  the  Center,  we  have  maintained  a  balanced  approach  of  physi¬ 
cal  modeling,  mathematical  analysis  and  numerical  computation.  The  applications  of  these  subjects  have 
served  as  an  important  guide  to  the  research  problems.  The  interaction  of  analysis  and  computation  has 
been  one  of  the  most  successful  features  of  the  Center. 

The  research  on  these  subjects  has  aimed  toward  the  following  goals: 

•  Development  of  mathematical  models 

•  Mathematical  analysis  of  physical  systems 

•  Development  and  analysis  of  numerical  methods 

•  Computation  of  particular  solutions 

•  Computation  of  solutions  when  analysis  is  impossible. 

These  goals  have  been  met  in  a  number  of  research  problems,  which  are  described  in  Sections  IV  and  VII. 
Two  key  elements  in  the  success  of  the  Center  have  been  the  strong  computational  environment  that  we 
have  assembled  and  the  high  quality  postdoctoral  visitors  and  graduate  students  that  we  were  able  to 
attract. 

Another  important  goal  has  been  to  generate  interest  in  these  research  subjects  among  students  and 
other  mathematicians  and  scientists.  This  has  been  accomplished  through  numerous  graduate  courses, 
seminars  and  publications  and  two  workshops. 


III.  Personnel 


The  core  of  the  Center  consists  of  4  faculty  members: 

Russel  Caflisch 
Graeme  Milton 
George  Papanicolaou 
Lu  Ting 

While  a  member  of  the  Center,  Graeme  Milton  was  awarded  the  prestigious  Packard  Fellowship  for  Sci¬ 
ence  and  Engineering,  as  well  as  a  Sloan  Fellowship.  Russel  Caflisch  received  a  Sloan  Fellowship. 

Postdoctoral  researchers  in  the  Center  have  included  the  following: 

Marco  Avellaneda  (1987-88) 

Thomas  Chacon  (1988) 

Tom  Hou  (1987-89) 

Robert  Krasny  (1987) 

Michael  Landman  (1987-89) 

John  Lowengrub  (1988-89) 

Robert  Palais  (1986-87) 

Demitri  Papageorgiou  (1986-87) 

Bruce  Pitman  (1986-87) 

Marie  Postel  (1988-89) 

Sophie  Weinryb  (1988-89) 

Marco  Avellaneda  joined  the  faculty  of  the  Courant  Institute  as  an  assistant  professor  in  Fall  1988.  Tom 
Hou  will  be  joining  the  faculty  of  the  Courant  Institute  as  an  assistant  professor  in  Fall  1989. 

Long-term  scientific  visitors  to  the  Center  have  included: 

James  Berryman  (Lawrence  Livermore  National  Laboratory) 

Alexander  Frenkel  (City  College,  New  York) 
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Christos  Flytzanis  (Ecole  Polytechnique,  France) 

Brenton  LeMesurier  (University  of  Arizona) 

Oscar  Orellana  (Santa  Maria  University,  Chile) 

Dale  Pullin  (University  of  Queensland,  Australia) 

Gregory  Sivashinsky  (City  College,  New  York) 

Pierre-Louis  Sulem  (University  of  Tel-Aviv,  Israel) 

Catherine  Sulem  (Ben-Gurion  University,  Israel) 

Graduate  students  (and  their  Ph.D.  date  and  research  project)  in  the  Center  have  included: 

Mahmoud  Affouf  (Ph.D.  1988,  computations  of  dynamic  phase  transitions) 

Mark  Asch  (pulse  reflection  from  random  media) 

Albert  Fannjiang  (nonlinear  random  media) 

Robert  Knapp  (Ph.D.  1988,  nonlinear  random  media) 

John  Lowengrub  (Ph.D.  1988,  vortex  methods) 

William  Morokoff  (numerical  methods  for  the  Boltzmann  equation) 

Andre  Nachbin  (Ph.D.  1989,  gravity  waves  in  shallow  water  with  irregular  bottom) 

Vincenzo  Nesi  (Ph.D.  1989,  optimal  bounds  for  composites) 

Oscar  Orellana  (Ph.D.  1987,  vortex  sheets) 

Gabriella  Puppo  (Ph.D.  1989,  computation  of  boundary  layers) 

Julia  Rennenkampf  (pulse  reflection  from  random  media) 

Michael  Siegel  (Ph.D.  1989,  singularities  in  Rayleigh-Taylor  flow) 

Xiao  Wang  (analysis  and  computation  of  singular  solutions  of  the  nonlinear  Schroedinger  equation) 
Xue  Xin  ( nonlinear  homogenization) 

Jing-Yi  Zhu  (Ph.D.  1989,  adaptive  vortex  method) 
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IV.  Synopses  of  Research  Results 

Here  we  give  synopses  of  several  outstanding  research  results  of  the  Center.  More  detailed  descrip¬ 
tions  of  these  and  other  research  projects  is  found  in  section  VII.  Some  of  the  most  important  results  of 
the  Center  have  been  the  following: 

1 .  Computations  and  numerical  analysis  of  the  vortex  method  for  vortex  sheets  were  carried  out  by 
Krasny  and  by  Caflisch  and  Lowengrub. 

2.  Exact  singular  solutions  of  the  vortex  sheet  equations  were  found  by  Caflisch  and  Orellana. 

3.  Caflisch  and  Siegel  have  developed  a  theory  for  the  formation  of  singularities  in  Rayleigh-Taylor 
flows,  which  gives  a  prediction  of  their  form  and  time  of  development. 

4.  The  stability  analysis  of  the  classical  three  point-vortex  problem  was  completed  by  Tavantzis  and 

Ting. 

5.  Ting  and  coworkers  have  analyzed  reflection  and  transmission  of  an  acoutic  wave  through  a  bub¬ 
bly  layer  in  a  liquid. 

6.  Convergence  of  the  point  vortex  method  for  inviscid,  incompressible  fluid  flow  was  proved  by 
Goodman,  Hou  and  Lowengrub.  This  new  result  is  expected  to  greatly  change  the  understanding  and  use 
of  vortex  methods. 

7.  An  adaptive  vortex  method  was  developed,  implemented  and  analyzed  by  by  Hou,  Papageorgiou, 
Papanicolaou  and  Zhu. 

8.  Papanicolaou  and  Zhu  have  developed  and  implemented  a  particle  method  for  computation  of 
dilute  fluid-particle  flows. 

9.  A  finite  element  vortex  method  was  developed,  implemented  and  analyzed  by  Chacon  and  Hou. 

10.  The  formulation  and  solution  of  inverse  problems  for  pulse  reflection  from  random  media  was 
developed  by  Papanicolaou  and  coworkers. 

11.  Nachbin  and  Papanicolaou  have  developed  and  implemented  a  computational  method  for  water 
waves  in  a  channel  with  an  irregular  bottom. 

12.  Kohn  and  Milton  have  found  generalizations  of  the  Hashin-Shtrikman  bounds  for  anisotropic 
materials. 

13.  Optimal  upper  and  lower  bounds  on  the  conductivity  of  a  polycrystalline  material  were  derived 


by  Milton  and  coworkers. 
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14.  Milton  has  derived  a  macroscopic  theory  of  the  Hall  effect  in  composite  materials. 

15.  A  numerical  homogenization  method  originally  developed  by  Engquist  and  Hou  for  hyperbolic 
equations  was  extended  by  Avellaneda,  Hou  and  Papanicolaou  to  elliptic  equations  in  several  dimensions. 

16.  The  form  of  the  focusing  singularity  of  the  nonlinear  Schroedinger  equation  was  determined  by 
Landman,  Le  Mesurier,  Papanicolaou,  Sulem  and  Sulem. 

17.  Localization  and  bistability  for  a  nonlinear,  random  optical  media  was  analyzed  and  simulated 
by  Knapp,  Papanicolaou  and  White. 
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V.  Computational  Environment 

The  computational  environment  has  been  a  very  important  component  of  the  Center.  With  funds 
from  the  AFOSR  and  the  NSF  and  with  substantial  contributions  from  New  York  University,  the  Center 
acquired  a  system  of  SUN  3  workstations  linked  by  an  ethemet,  a  Celerity  mini  computer,  an  Iris  graphics 
workstation  and  an  array  of  peripherals  and  software. 

Ease  of  use  of  the  workstations  and  general  availability  of  computational  cycles  has  transformed  the 
role  of  computations  in  our  research  efforts.  In  addition  to  their  traditional  use  in  providing  highly  accu¬ 
rate  solution  of  well-formulated  systems,  computations  now  serve  as  an  experimental  mathematical  tool 
for  investigating  phenomena  that  are  poorly  understood  and  for  testing  hypothetical  properties.  This 
experimental  role  for  computations  has  been  particularly  strong,  for  example,  in  the  our  investigations  of 
singularity  formation  for  the  nonlinear  Schroedinger  equation,  singularity  lormation  for  vortex  sheets, 
and  the  solution  of  inverse  problems  in  random  media. 

The  principal  hardware  components  of  our  computer  system  are 

1  Sun-3/180  fileserver  and  color  workstation  (the  main  server  at  present) 

10  Sun-3/50  and  1  SUN-3/60  (purchased  with  NSF  funds)  workstations 

1  SUN-3/280  server  (purchased  with  CIMS  funds  but  not  fully  equipped  at  present) 

1  Celerity  1200  mini  computer  (purchased  with  NSF  funds  in  1985) 

1  Iris  4D/50-GT  graphics  workstation. 

In  addition  to  the  standard  UNIX-based  software  distributed  with  the  machines,  we  have  several 
public-domain  software  packages  (cmlib,  auto,  near  graphics,  etc.)  and  the  following  software  items  that 
we  purchased: 

Precision  Visuals  DI3000,  PICSURE  and  Metafile  translator  (graphics  package) 

NAG  library  of  scientific  routines. 

Graphics  animation  system  (GAS)  and  Surface  plot  system  (SURF)  for  the  IRIS,  from  Sterling 
Software 

Our  supercomputing  is  done  at  the  John  von  Neumann  Center  at  Princeton  and  is  sponsored  by  the 
NSF.  We  have  a  1.5  megabit/second  communications  link  to  Princeton’s  Cyber  205’s  and  ETA  10’s.  The 
speed  of  data  transmission  to  and  from  Princeton  make  it  the  preferred  supercomputing  center  for  us.  We 
have  accounts  at  San  Diego,  NASA-Ames,  Kirtland  AFB  and  Pittsburgh,  but  the  links  are  slow  and  so  we 
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do  not  use  them  very  much,  even  though  they  have  Cray’s  running  UNIX  making  them  much  more  con¬ 
venient  for  us. 


fi 
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VI.  Workshops 

Two  workshops  were  organized  by  the  Center: 

Workshop  on  Mathematical  Aspects  of  Vortex  Dynamics 

The  Workshop  on  Mathematical  Aspects  of  Vortex  Dynamics  was  held  at  the  Xerox  Training 
Center  in  Leesburg,  Virginia,  on  April  25  -  27,  1988.  The  Center  and  the  Society  for  Industrial  and 
Applied  Mathematics  (SIAM)  were  joint  sponsors.  The  workshop  brought  together  researchers  in 
mathematical  analysis,  scientific  computing  and  theoretical  physics  and  engineering,  from  universities, 
industry  and  government  research  laboratories,  who  are  working  on  vortex  dynamics  in  fluids.  A  large 
variety  of  fluid  problems  were  discussed,  such  as  vortex  sheets  and  shear  layers,  the  vortex  reconnection 
problem,  vortex  methods,  coherent  vortex  structures,  fluid  dynamic  stability,  systems  of  point  vortices, 
turbulence  modeling  and  other  applications.  The  results  presented  included  mathematical  existence 
theorems,  results  on  convergence  of  numerical  methods,  computations  and  qualitative  analysis  of  flows. 

Topics  and  speakers  for  the  meeting  were  selected  by  a  Program  Committee  consisting  of  Hassan 
Aref  (UC  San  Diego),  Gregory  Baker  (Exxon  Corporate  Research  and  Ohio  State  University),  J.  Thomas 
Beale  (Duke  University),  Russel  Caflisch  (chairman,  Courant  Institute),  Steven  Childress  (Courant  Insti¬ 
tute),  and  Sidney  Lcibovich  (Cornell  University). 

There  were  approximately  60  attendees  at  the  workshop,  representing  universities,  government  and 
industry.  The  speakers  and  the  titles  of  their  presentations  were  the  following: 

1.  C.  Anderson,  (University  of  California,  Los  Angeles),  "An  Investigation  of  the  Accuracy  of  Vorti- 
city  Boundary  Condition" 

2.  B.  J.  Bayly,  (Courant  Institute),  "Computations  of  Broad-Band  Instabilities  in  a  Class  of  Closed- 
Streamline  Flows" 

3.  M.  Berger  (U.  Massachusetts),  "Remarks  on  Vortex  Breakdown" 

4.  T.  Buttke,  (Princeton  University),  "Turbulence  in  Super  Fluids" 

5.  L.  J.  Campbell,  (Los  Alamos  National  Laboratory),  "Vortex  Lattices  in  Theory  and  Practice" 

6.  G.H.  Cottet,  (UCLA  and  University  of  Paris),  "New  Boundary  Conditions  for  Navier  Stokes  Equa¬ 
tions  in  the  Vorticity  Formulation" 

7.  D.  Dritschel,  (Cambridge  University),  "Stripping  of  Two-Dimensional  Vortices  by  Weak  Strain  and 
Shear" 

8.  T.  Kambe,  (University  of  Tokyo),  "Analysis  of  Vortex-Sheet  Motions  Based  on  the  Theory  of 
Hyperfunctions" 

9.  A.  Leonard  and  G.  Winckelmans,  (California  Institute  of  Technology),  "Improved  Vortex  Methods 
for  Three-Dimensional  Flows  with  Application  to  the  Interactions  of  Two  Vortex  Rings" 

10.  D.G.  Long,  (Princeton  University),  "The  Random  Vortex  Method  in  Two  Dimensions" 
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11.  J.  Lowengrub,  (Courant  Institute),  "Convergence  of  the  Vortex  Method  for  Vortex  Sheets" 

12.  J.  L.  Lumley,  (Cornell  University),  "The  Dynamics  of  Vortex  Structures  in  the  Wall  Region  of  a 
Turbulent  Boundary  Layer" 

13.  T.  S.  Lundgren,  (University  of  Minnesota),  "A  Free  Surface  Vortex  Method  With  Small  Viscous 
Effects" 

14.  F.  S.  Marcus,  (University  of  California,  Berkeley),  "Vortex  Dynamics  in  Shearing  Flows" 

15.  D.  I.  Meiron,  (California  Institute  of  Technology),  S.  A.  Orszag,  M.  J.  Shelley,  (Princeton  Univer¬ 
sity),  "A  Numerical  Study  of  Vortex  Reconnection" 

16.  Richard  Peltz,  (Rutgers  University),  "Spectral  Methods  on  Hypercube  Computers" 

17.  C.  Pozrikidis,  (University  of  California  at  San  Diego),  "Two-Dimensional  Vortex  Dynamics  in 
Flows  of  Mixed  Shear  Layer-Wake  Type" 

18.  D.  I.  Pullin,  (The  University  of  Queensland,  Australia),  "On  Similarity  Solutions  for  the  Self- 
Induced  Motion  of  Two-Branched  Vortex  Sheets" 

19.  M.  Pulvirenti,  (Universita  dell’Aquila,  Italy),  "On  the  Invariant  Measures  for  the  2-D  Euler  Flow" 

20.  M.  J.  Shelley  and  C.  R.  Baker,  (Princeton  University),  "On  the  Connection  between  Thin  Vortex 
Layers  and  Vortex  Sheets" 

21.  G.  Tryggvason,  (The  University  of  Michigan),  "Vortex  Dynamics  of  Stratified  Flows" 

22.  B.  Turkington,  A.  Eydeland,  (University  of  Massachusetts),  "An  Iterative  Method  For  Computing 
Steady  Vortex  Flow  Systems" 

The  proceedings  of  the  workshop  will  be  published  in  the  SIAM  Conference  Proceedings  Series. 


Workshop  on  Random  Media  and  Composites 

Composite  and  random  media  include  a  wide  variety  of  materials  ranging  from  alloys,  emulsions 
and  colloidal  suspensions  to  polycrystalline  aggregates,  porous  rocks  and  fiberous  substances.  The  Air 
Force  URI  Center  at  the  Courant  Institute  and  the  Society  for  Industrial  and  Applied  Mathematics 
(SIAM)  held  a  conference  on  random  media  and  composites  from  December  7th  to  10th  1988  inclusive, 
at  the  Xerox  Training  Center  in  Leesburg,  Virginia.  The  main  emphasis  of  the  conference  was  on  theoret¬ 
ical  developments  in  understanding  the  effective  moduli  of  composites  and  on  localization  in  random 
media.  The  scope  also  extended  to  percolation,  to  non-linear  effects  and  micromechanics  and  to  compo¬ 
site  effects  in  superconductivity. 

The  organizing  committee  for  the  scientific  program  of  the  conference  consisted  of  Graeme  Milton 
(Courant  Institute),  chairman;  David  Johnson  (Schlumberger-Doll  Research  Laboratories);  Robert  Kohn 
(Courant  Institute);  George  Papanicolaou  (Courant  Institute);  Pabitra  Sen  (Schlumberger-Doll  Research 
Laboratories);  Ping  Sheng  (Exxon  Corporate  Research);  Luc  Tartar  (Camegie-Mellon  University);  Benja¬ 
min  White  (Exxon  Corporate  Research).  The  committee  represented  a  balance  of  research  scientists  from 
universities  and  industry.  Johnson,  Milton,  Sen  and  Sheng  have  close  connections  with  the  physics  com¬ 
munity  while  Kohn,  Papanicolaou,  Tartar,  and  White  are  most  familiar  with  recent  mathematical  develop¬ 
ments. 

There  were  approximately  80  attendees  at  the  workshop,  representing  universities,  government  and 
industry.  The  speakers  and  the  titles  of  their  presentations  were  the  following: 

1.  A.  Acrivos,  (City  University  of  New  York,  New  York),  "An  Effective  Continuum  Theory  for 
Estimating  Transport  Parameters  in  Two-Phase  Random  Media" 

2.  M.  Avellaneda,  (Courant  Institute),  "Optimal  Bounds  and  Microgeometries  for  Polycrystalline 
Aggregates" 

3.  R.  Barrera,  (Universidad  Nacional  de  Mexico),  "An  Intuitive  Approach  and  a  Diagrammatic  Formu¬ 
lation  of  the  Renormalized  Polarizability  in  the  Maxwell-Gamett  Theory" 

4.  D.  Bergman,  (Tel-Aviv  University,  Israel),  "Nonlinear  Behavior  and  1/f  Noise  Near  a  Conductivity 
Threshold-The  Importance  or  Local  Geometry  and  the  Failure  of  the  Effective  Medium  Theory" 

5.  P.M.  Duxbury,  (Michigan  State  University),  "Scaling  Theory  and  Fluctuations  in  the  Strength  of 
Composites" 

6.  D.  Fisher,  (Princeton  University),  "Onset  of  Superfluidity  in  Random  Media" 

7.  K.  Golden,  (Princeton  University),  "Convexity  in  Random  Resistor  Networks" 

8.  L.  Greene,  (Bell  Communications  Research),  "High  Temperature  Superconductivity" 

9.  B.  Halperin,  (Harvard  University),  "Percolation  and  Transport  in  Highly  Disordered  Systems" 

10.  R.  James,  (University  of  Minnesota,  Minneapolis),  "Theory  of  Martensitic  Transformations" 

11.  A.  Lagendijk,  (University  of  Amsterdam,  The  Netherlands),  "Experiments  on  the  Localization  of 
Light" 


12.  C.  Lobb,  (Harvard  University),  "Simulations  of  Composites:  From  Ohmic  Response  to  the  Break¬ 
down  of  Superconductivity" 

13.  J.  Machta,  (University  of  Massachusetts),  "Localization  of  the  Sound  Modes  of  Helium  in  Porous 
Media  and  Symmetry  Breaking  by  a  Row" 

14.  J.  Maynard,  (Pennsylvania  State  University),  "Acoustic  Anderson  Localization" 

15.  W.  Means,  (State  University  of  New  York),  "Synkinematic  Microscopy  of  Transparent  Polycrys¬ 
tals" 

16.  F.  Murat,  (University  of  Paris  VI,  France),  "Homogenization  in  Thermoelasticity" 

17.  T.  W.  Noh,  (Cornell  University),  "Infrared  Electrodynamic  Properties  of  High  Tc  Superconductors" 

18.  L.  Poladian,  (University  of  Sidney,  Australia),  "Critical  Behaviour  of  Effective  Moduli  Determined 
from  Nearest  Neighbour  Interactions” 

19.  S.  Redner,  (Boston  University),  "Transport  and  Dispersion  in  Random  Media" 

20.  J.  Rubinstein,  (Israel  Institute  of  Technology),  "Bounds  on  Various  Electrostatic  and  Hydrodynamic 
Capacities" 

21.  D.  Stroud,  (Ohio  State  University),  "Theory  of  the  Transport  and  A.C.  Properties  of  Superconduct¬ 
ing  Composites" 

22.  W.  Symes,  (Rice  University),  "Propagation  of  Regularity  and  Inverse  Problems" 

23.  D.  Tanner,  (University  of  Florida,  Gainsville),  "Infrared  Studies  of  High-Tc  Materials" 

24.  B.  White,  (Exxon  Research  and  Engineering  Co.),  "Pulse  Reflection  from  a  Randomly-Stratified 
Medium" 

25.  John  Willis,  (University  of  Bath,  England),  "The  Overall  Behaviour  of  a  Nonlinear  Matrix  Rein¬ 
forced  by  Rigid  Inclusions  or  Weakened  by  Cavities" 

26.  P.Z.  Wong,  (University  of  Massachusetts),  "Transport  Properties  of  Real  and  Artificial  Rocks" 
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VII.  Descriptions  of  Research  Projects 


VH.l  Fluid  Dynamics 

(1)  Analysis  and  Computation  for  Vortex  Sheets  (Caflisch,  Krasny,  Lowengrub,  Orellana) 

A  vortex  sheet  in  two-dimensional,  inviscid,  incompressible  fluid  flow  is  a  curve  along  which  the 
velocity  is  discontinuous  in  its  tangential  component.  Vortex  sheets  are  generated  in  a  variety  of  labora¬ 
tory  and  natural  flows,  such  as  mixing  flow  past  a  splitter  plate,  separated  boundary  layers,  and  flow  over 
delta  wings.  The  instability  of  vortex  sheets  is  one  of  the  primary  mechanisms  for  generation  of  small 
scales  and  the  onset  of  turbulence  in  high  Reynolds  number  flows.  In  fact  the  growth  rate  for  this  insta¬ 
bility  can  be  arbitrarily  large  since  it  is  proportional  to  the  wavenumber. 

This  research  project  is  on  computation  and  analysis  of  this  strongly  unstable  (even  ill-posed)  prob¬ 
lem.  In  particular  the  vortex  sheet  develops  a  singularity  (a  point  at  which  the  curvature  of  the  sheet  is 
infinite)  and  then  rolls-up  into  a  spiral.  The  main  research  goals  are  the  following: 

(1)  Determination  of  a  mathematical  formulation  in  which  the  vortex  sheet  problem  is  well-posed. 

(2)  Development  of  numerical  methods  for  the  vortex  sheet  problem,  analysis  of  convergence  of  these 
methods,  and  investigation  of  particular  physical  phenomena  by  numerical  computation. 

(3)  Construction  of  examples  of  singularity  formation  on  vortex  sheets. 

(4)  Determination  of  the  generic  form  for  singularity  formation. 

(5)  Investigation  of  the  behavior  of  the  vortex  sheet  after  the  singularity  formation,  including  roll-up. 

As  described  below,  (1-3)  are  fairly  well  completed;  while  we  have  only  partial  results  on  (4)  and  (5). 
This  project  has  been  initiated  by  the  numerical  investigations  of  Krasny[8,9],  All  of  the  analysis  has 
been  inspired  by  his  computational  results. 

One  of  the  main  mathematical  motivations  for  this  work  is  the  importance  of  singularities  for  a 
variety  of  fluid  problems.  Most  important  is  the  possible  development  of  singularities  in  three- 
dimensional  inviscid  flow,  which  would  be  related  to  the  development  of  turbulence.  On  the  other  hand, 
in  two  dimensions  singularities  do  not  develop  from  smooth  initial  data.  The  vortex  sheet  problem  is  of 
intermediate  difficulty:  singular  initial  data  (a  vortex  sheet)  in  two  dimensions  are  found  to  become  more 
singular  after  some  time.  An  even  stronger  singularity,  involving  energy  concentration,  has  been 
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proposed  for  two  dimensional  flows  by  DiPema  and  Majda  [4-6], 

The  vortex  sheet  is  described  by  a  complex  function  z  (y,t)  =  x  (y,t)  +  iy  (y,t)  in  which  y  is  the  circu¬ 
lation  variable  and  z  is  the  complex  position  of  the  sheet.  The  evolution  of  the  sheet  is  governed  by  the 
Birkhoff-Rott  equation 

oo 

a,*"M  =  B[z]  =  (2ni)-'  PV  j  (z  (y,t)  -  Z  (Y  '.t))"1  dy '  (1) 

in  which  the  integral  is  a  Cauchy  principal  value  integral.  An  exact  steady  solution  of  this  equation  is 
z  =  Y  corresponding  to  a  flat  vortex  sheet  of  uniform  strength.  For  a  nearly  flat  vortex  sheet  the  solution 
can  be  written  as  z  -y+s(y,t)  with  s  small.  The  principal  linear  modes  for  s  are 
s-t*  =  (1±0  exp{iky±kt!2).  This  Kelvin-Helmholtz  instability  has  a  linear  growth  rate  that  is  propor¬ 
tional  to  the  wavenumber  k  and  hence  can  be  arbitrarily  large. 


Numerical  Computations  of  Vortex  Sheets 

Krasny  [8,9]  developed  two  methods  for  computing  the  evolution  of  vortex  sheets.  The  main 
numerical  difficulty  of  this  problem  is  that  the  physical  instability  of  the  vortex  sheet  can  magnify  the 
roundoff  errors  of  the  computation. 

In  the  first  method  the  roundoff  errors  are  controlled  by  spectral  filtering  at  a  fixed  filter  level  ey.  At 
each  step  of  the  computation  the  fourier  component  s(k)  is  set  to  zero  if  its  magnitude  is  below  ef.  In  this 
way  the  high  wavenumber  components  grow  only  if  the  nonlinearity  of  the  problem  generates  them.  This 
method  works  up  to  the  time  of  singularity  formation  and  gives  predictions  on  the  time  and  type  of  singu¬ 
larity  formation.  Meiron,  Baker  and  Orszag  [10]  also  computed  the  form  and  time  of  the  singularity  by 
an  alternative  method. 

Krasny’s  second  method  is  a  desingularization  of  the  Birkhoff-Rott  equation  (1).  Replace  the 
integral  operator  B  [z  ]  in  (1)  by 


S6[z](y)  =  —  PV  j - z(7>~2(Y') —  dy' 

2 ni  jL  I  z  (y>-z  (Y  0 1 2  +  82 


(2) 


which  equals  B  [z  ]  for  8  =  0  and  has  no  singularities  for  8  >  0.  Using  this  equation  Krasny  was  able  to 
compute  past  the  singularity  time  and  find  rolled-up  spiral  solutions,  as  well  as  interactions  of  these 
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spirals. 

Well-Posedness  of  Vortex  Sheets  for  Analytic  Initial  Data 

Kelvin-Helmholtz  instability  makes  the  vortex  sheet  problem  ill-posed  in  a  Sobolev  function  space, 
but  the  problem  is  well-poser  in  the  space  of  analytic  functions.  The  reason  for  this  is  that  the  linearized 
Birkhoff-Rott  equation  is  hyperbolic  in  the  imaginary  y  direction.  Thus  details  of  the  solution  move  on 
characteristics  at  speed  ±i/2.  Using  this  description  of  the  problem,  Caflisch  and  Orellana  [1]  proved  the 
following  long  time  existence  theorem: 

Suppose  the  initial  data  s  (r=0,  Y)  is  analytic  and  small  in  a  strip  {\Imy\  <  p).  Then  the  Birkhoff- 
Rott  equation  (1)  has  an  analytic  solution  z  =  y+  s(t,y)  for  a  time  interval  0  <  t  <  (2-a)p,  in  which  a  is  a 
small  nonzero  number. 

This  result  is  nearly  optimal,  since  it  is  valid  almost  up  to  the  observed  singularity  time.  An  earlier 
short-time  version  of  this  result  for  somewhat  more  general  initial  data  was  proved  by  Sulem,  Sulem, 
Bardos  and  Frischf  13].  The  proof  is  based  on  the  abstract  Cauchy-Kowalewski  theorem. 

Next  consider  Krasny’s  desingularized  equation  with  time  discretization  size  At  and  with  simulated 
roundoff  error  er,  given  by 

z(  (n  +l)Ar)  =  z(nAr)  +  Ar  B^zf/iAr)]  +  er,  (3) 

in  which  is  a  discrete  sum  approximation  to  the  integral  J5g.  Using  the  same  methods,  Caflisch  and 
Lowengrub[3]  proved  the  following  convergence  theorem  for  Krasny’s  numerical  method: 

Suppose  as  above  that  the  initial  data  s(r= 0,  y)  is  analytic  and  small  in  a  strip  {\Imy\  <  p}.  Also 
assume  that  the  roundoff  error  satisfies 


\er  I  <  At  e~1/Ax.  (4) 

Then  I  z  -  z  I  <  c  (At  +  Ax  +  8  +  e). 

This  theorem  is  only  for  a  short  time,  much  shorter  than  the  singularity  time.  The  restriction  (4)  on 


the  roundoff  error  size  is  the  numerical  analysis  interpretation  of  the  restriction  to  analytic  functions. 
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Singularities  •  Examples  and  the  Generic  Form  of  Singularities 

Singular  solutions  of  the  Birkhoff-Rott  equation  (1)  were  constructed  by  Caflisch  and  Orellana  [2] 
and  by  Duchon  and  Robert  [7].  These  solutions  are  analytic  in  y  for  t  =  0,  and  at  the  critical  time  t  =  tc 
they  have  the  form  z  =  y +  e  sgn{ y)  I  yl a  for  y  =  0  in  which  e  is  small.  Such  solutions  can  be  produced  for 
any  a  >  1;  i.e.  singularities  can  be  produced  in  anything  stronger  than  the  first  derivative.  In  particular 
the  curvature  of  the  sheet  will  be  infinite  at  the  singularity. 

These  solutions  are  not  believed  to  be  typical.  An  asymptotic  analysis  by  Moore  [11,12]  suggests 
that  the  typical  form  of  singularities  has  a  =  3/2.  The  reason  is  that  the  singularity  occurs  as  an  envelope 
of  characteristics,  and  on  such  an  envelope  the  typical  behavior  is  square  root  type.  A  precise  mathemati¬ 
cal  analysis  of  this  question  is  currently  being  carried  out  by  Caflisch,  Ercolani,  Hou  and  Semmes.  Their 
analysis  is  based  on  an  approximation  of  the  Birkhoff-Rott  equation  by  a  system  of  nonlinear  hyperbolic 
(in  the  imaginary  y  direction)  equations  of  the  form 

d,F  =M(F)  dyF.  (5) 

The  vector  F  is  ( s ,  s*,  Hs,  Hs*)  in  which  s(y)-s  (yi)  and  Hs  is  the  Hilbert  transform  of  s. 

The  vortex  sheet  singularity  and  subsequent  roll-up  of  the  sheet  act  as  a  strong  source  of  noise  in  a 
slightly  compressible  fluid.  An  analysis  of  the  resulting  sound  waves  is  found  in  [14], 

Caflisch  is  currently  writing  a  monograph  on  the  mathematical  and  numerical  aspects  of  vortex 
dynamics. 
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(2)  Singularity  Formation  for  Rayleigh-Taylor  Flows  (Caflisch,  Siegel) 

The  Rayleigh-Taylor  problem  consists  of  two  fluids  of  different  densities  pi  and  P2  separated  by  an 
interface  and  acted  on  by  gravity  g  or  some  other  external  force.  This  is  a  fundamental  problem  of  fluid 
dynamics  that  arises  in  a  large  variety  of  applications,  such  as  a  blast  wave  propagating  through  an  inter¬ 
face. 

A  steady  solution  of  the  problem  is  for  the  interface  to  be  flat  and  stationary  with  a  hydrostatic  pres¬ 
sure  distribution.  However  this  solution  is  unstable,  if  the  top  density  pi  is  larger  than  the  bottom  density 
p2.  As  a  result  fingers  of  the  heavier  top  fluid  are  pushed  into  the  lighter  bottom  fluid.  The  resulting  type 
of  fingering  is  strongly  dependent  on  the  value  of  the  Atwood  number  A,  defined  by  A  =  (pi-p2)/pi 

As  the  finger  pushes  into  the  bottom  fluid,  vorticity  is  baroclinically  generated  on  the  interface  so 
that  the  interface  is  unstable  due  to  Kelvin-Helmholtz  instability.  Thus  if  surface  tension  and  viscosity 
are  ignored,  we  expect  singularities  to  form  along  the  interface,  followed  by  roll-up  of  the  interface,  as 
described  in  subsection  (1).  Vorticity  is  of  opposite  sign  on  the  two  sides  of  a  finger,  so  that  singularities 
of  opposite  sign  are  expected  to  form  on  the  two  sides.  This  process  is  important  as  a  mechanism  for  gen¬ 
eration  of  small  scale  structures  in  the  flow.  If  surface  tension  or  viscosity  were  small  but  nonzero,  the 
singularity  would  be  slightly  smoothed-out,  but  the  generation  of  small  scales  would  still  be  present. 

This  roll-up  has  been  experimentally  observed  and  numerically  computed  [1].  Singularities  are 
more  difficult  to  compute.  But  in  the  Boussinesq  approximation  with  A— >0,  Ag— >1,  singularities  have 
been  computed  by  Pugh  [5],  a  student  of  Derek  Moore.  At  the  other  extreme,  A  =  1,  it  is  observed  that 
the  positions  of  the  singularities  on  the  two  sides  of  the  finger  tend  to  merge  at  the  tip.  Since  the  two 
singularities  come  with  opposite  sign,  it  has  been  conjectured  [1]  that  they  will  cancel  and  that  the  flow 
for  A  =  1  will  be  smooth. 

In  his  Ph.D.  thesis  under  the  direction  of  Caflisch,  Siegel  [6]  has  developed  a  simple  approximate 
theory  to  describe  the  Rayleigh-Taylor  problem.  This  theory  is  a  particular  case  of  a  more  general  approx¬ 
imation  method  developed  in  [2].  Whereas  the  Rayleigh-Taylor  equations  are  nonlocal  integro- 
differential  equations,  the  approximate  theory  consists  of  a  system  of  first  order  pde’s.  For  a  nearly  flat 
interface  these  pde’s  can  be  approximately  solved  by  perturbation  theory.  A  singularity  is  seen  to  form  as 
an  envelope  of  the  characteristics  of  the  pde,  analogous  to  the  singularity  for  the  Kelvin-Helmholtz  equa¬ 
tion  described  by  Moore  [3,4].  In  this  way  Siegel  finds  an  expression  for  the  time  and  type  of  singularity 
formation,  the  first  such  analytic  expression  for  this  problem.  This  result  completely  agrees  with  Pugh’s 
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numerical  results  for  the  Boussinesq  limit. 

Siegel  is  currently  analyzing  the  cancellation  of  the  singularity  for  A  =  1  and  carrying  out  more 
extensive  computations  to  verify  his  analytic  results.  It  may  also  be  possible  to  include  surface  tension 
effects. 
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(3)  Dynamics  of  Point  Vortices  (Ting) 

The  motion  of  N  free  point  vortices  in  a  two  dimensional  incompressible  inviscid  fluid  and  their 
induced  flow  field  has  been  a  subject  of  interest  in  studies  of  Hamiltonian  dynamical  systems,  numerical 
simulations  of  chaotic  motions  and  physical  applications  (see  e.g.,  the  recent  review  paper  by  Aref  [1] 
and  the  references  therein).  The  dynamics  of  three  point  vortices  was  studied  by  Groebli  [2]  (1877)  and 
by  Synge  [3]  (1949).  Unaware  of  the  earlier  works,  Novikov  [4]  (1975)  and  Aref  [5]  (1979)  rediscovered 
the  results  obtained  by  Groebli  and  Synge.  In  Synge’s  analysis  [3],  the  lengths  of  the  sides  of  the  triangle 
formed  by  the  three  vortices  were  used  as  prime  variables.  The  critical  states  at  which  the  lengths  of  the 
sides  remain  fixed  throughout  the  motion  were  found  to  be  either  equilateral  triangles  or  collinear 
configurations.  The  equilateral  configurations  were  shown  to  be  stable  or  unstable  depending  on  whether 
the  sum  of  the  products  of  strengths  K  was  greater  or  less  than  zero  respectively.  In  the  case  K  =  0,  a  one 
parameter  family  of  solutions  of  contracting  and  another  of  expanding  similar  triangles  were  found. 

Ting  and  Tavantzis  [6]  continued  Synge’s  analysis  so  as  to  obtain  the  global  behavior  of  the  trajec¬ 
tories  representing  the  configurations  formed  by  the  three  vortices.  They  are  able  to  completely  classify 
the  trajectories  given  the  strengths  of  the  three  vortices  and  describe  the  trajectory  qualitatively  once  the 
initial  configuration  is  specified.  They  show  that  for  the  special  case  of  K  =  0,  the  family  of  contracting 
similar  solutions  is  always  unstable  while  the  family  of  expanding  ones  is  stable.  They  then  study  the  crit¬ 
ical  states  for  collinear  configurations  in  the  general  case  where  K  is  greater  or  less  than  0,  and  show  that 
there  are  either  six  or  four  critical  states  depending  on  the  strengths  of  the  vortices.  When  there  arc  six 
collinear  critical  states,  three  of  them  are  always  stable,  one  is  not,  while  the  remaining  two  are  unstable 
(stable)  if  the  equilateral  triangle  configuration  is  stable  (unstable).  When  there  are  only  four  collinear 
critical  states,  they  all  are  stable  while  the  equilateral  triangle  configuration  is  always  unstable.  Since 
there  are  two  equilateral  triangle  configurations,  clockwise  and  counter-clockwise  arrangements  of  the 
three  vortices,  the  sum  of  the  indices  of  all  the  critical  states  is  equal  to  +2  regardless  of  whether  K  is 
greater  or  less  than  0.  An  integral  invariant  in  trilinear  coordinates  is  derived.  Knowing  all  the  critical 
points  and  the  integral  invariant,  we  obtain  the  global  behavior  of  the  trajectories. 

References  are  listed  at  the  end  of  subsection  (6). 
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(4)  Motion  and  Interaction  of  Vortices  with  Diffusive  Core  Structures  (Ting) 

In  April  25-29,  1988  ten  invited  lectures  were  given  by  L.  Ting  at  the  Acrodynamisches  Institut, 
Rheinisch-Westfaelische  Technische  Hochschulc  Aachen  in  commemoration  of  the  75th  anniversary  of 
the  founding  of  the  Institute.  The  lecture  series,  entitled  "Motion  and  Interaction  of  Vortices  with  Dif¬ 
fusive  Core  Structures",  is  based  on  the  past  research  of  Ting  and  his  coworicers  sponsored  by  AFOSR, 
ONR  and  NASA  since  1964  (see  references  7-17).  The  outlines  of  the  ten  lectures  are  : 

(1)  The  length  and  time  scales  of  a  vortex  dominated  flow,  governing  equations,  the  integral  invariants, 
the  far  field  behavior  and  the  approximate  boundary  conditions  for  a  computational  domain  [7], 
(Lecture  I,  II) 

(2)  The  classical  inviscid  theory,  motion  of  a  spinning  disc,  the  Biot-Savart  formula  and  the  singular 
terms;  diffusion  of  an  isolated  vortex  and  the  optimum  similarity  solution  [8]. (Lecture  III) 

(3)  Matched  asymptotic  solutions  of  two  dimensional  vortical  flows,  multiple  length  and  time  scales, 
velocity  of  a  vortex  center  and  its  diffusive  core  structures,  the  optimum  similarity  solutions;  the 
meaning  of  the  leading  and  higher  order  solutions  [9].  (  Lecture  IV,  V) 

(4)  Motion  of  slender  vortex  filaments,  the  coupling  between  the  axial  and  circumferential  velocity 
components  in  the  core  structure,  the  physical  meaning,  the  ma'ched  asymptotic  solutions,  interac¬ 
tion  of  vortex  filaments  [10-12],  (Lecture  VI,  VII) 

(5)  Practical  limits  of  the  asymptotic  solutions,  studies  of  the  merging  of  vortices,  the  rules  of  merging 
[13-15],  (Lecture  VIII) 

(6)  Approximate  solutions  of  N-S  equations  using  vortices  with  diffusive  core  structures  [16,  17]. 
(l.ecture  IX) 

(7)  Concluding  remarks,  discussions.  (Lecture  X) 

An  extended  version  of  these  lectures  is  being  prepared  for  publication  as  a  monograph  in  colla¬ 
boration  with  Rupert  Klein,  who  was  at  Aachen  and  is  now  a  visiting  scholar  in  the  Department  of 
Mathematics,  Princeton  University. 
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(5)  Noise  Generation  from  Vortical  Flows  (Ting) 

A  primary  acoustic  source  in  unsteady  flows  is  the  time-dependence  of  the  vorticity.  Analysis  of 
this  source  of  noise  generation  is  of  importance  in  a  wide  variety  of  applications,  such  as  airplane  noise. 

Far  field  noise  induced  by  unsteady  vortical  flows  has  been  studied  by  many  investigators  interested 
in  sound  generation  by  turbulence.  In  his  pioneering  work,  Lighthill  [18]  found  that  the  leading  order 
acoustic  field  should  be  quadrupoles.  Later  Ribner  [19]  pointed  out  that  there  could  be  an  acoustic  source 
due  to  fluid  dilatation.  An  asymptotic  analysis  has  been  done  by  Crow  [20]  where  he  matched  the  inner 
low  Mach  number  inviscid  flow  to  an  outer  acoustic  field  and  showed  that  for  an  inviscid  isentropic  flow, 
the  leading  order  acoustic  field  could  be  composed  of  quadrupoles  only. 

In  the  investigation  of  Miksis  and  Ting  [21]  we  study  the  acoustic  field  induced  by  an  unsteady 
vortical  flow  with  viscosity.  We  consider  the  limit  where  the  reference  length  scale  of  the  vortical  flow  is 
much  smaller  than  the  reference  length  scale  of  the  acoustic  field.  At  low  Mach  number,  Af,  the  leading 
order  solution  of  the  vortical  flow  is  an  incompressible  viscous  flow.  We  match  the  far  field  behavior  of 
the  flow  [7]  to  the  quadrupoles  of  the  acoustic  field  and  relate  their  strengths  to  linear  combinations  of 
second  moments  of  vorticity.  The  compressibility  effect  of  the  vortical  flow,  of  0{M2),  induces  an 
acoustic  source  which  is  of  the  same  order  as  the  quadrupoles.  We  also  relate  the  strength  of  the  acoustic 
source  to  the  total  dissipation  energy  of  the  incompressible  vortical  flow. 

The  leading  order  acoustic  pressure  is  0(M3)  and  is  composed  of  acoustic  quadrapoles  and  an 
acoutic  source,  i.e. 
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where  x  is  the  retarded  time  t-r.  The  first  term  represents  the  quadropoles  with  strength 

B  (6,  y.x)  =  (9Ft(x)ijik  +  2^(Hj(x)  -  Hk(  x))). 


with  i.jjc  in  cyclic  order.  Here  £,,  i=l ,2,3  denote  the  three  components  of  x/r  and  arc  defined  by  the 
spherical  coordinates  0  and  y.  The  strengths  F,,  H,  are  related  to  the  second  moments  of  vorticity  in  the 


inner  region,  i.e. 

F,(0  =  <co i(x)  -  x\)>  and  //, (0  =  <2co tXjXk  -  (oyx*x,  -  0)*x;x; > 


for  i=l,2,3  with  i  j,k  in  cyclic  order.  We  note  that  H  j  +  H2+H 3  =  0,  so  that  there  are  only  five  quadro- 
poles. 


The  second  term  in  the  acoustic  pressure  (1)  represents  an  acoustic  source  with  strength  m( t) 
related  to  the  total  dissipation  energy  of  the  inner  region,  i.e. 
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These  results  are  applied  to  turbulent  flows  and  their  implications  for  the  mean  flow  are  discussed. 

It  is  known  that  the  far  field  behavior  of  the  vector  velocity  potential  induced  by  an  exponentially 
decaying  vorticity  distribution  is  given  by  a  power  series  in  r-1  where  r  is  the  distance  from  the  origin 
[7],  Now  we  can  show  that  the  nth  term  in  the  power  series  for  the  far  field  velocity  depends  only  on 
2n+l  linear  combinations  of  the  nth  moments  of  vorticity,  for  which  there  are  (3/2)(n+2)(n+l)  of  them. 
This  investigation  is  being  prepared  for  publication  by  R.  Klein  and  L.  Ting  [22]. 


References  are  listed  at  the  end  of  subsection  (6). 


-27- 


(6)  Wave  Propagation  and  Multiphase  Flow  (Ting) 

In  continuing  our  past  research  on  wave  propagation  through  bubbly  media,  we  investigated  the 
transmission  and  reflection  of  an  acoustic  wave  through  a  bubbly  layer.  We  use  the  nonlinear  model 
equations  of  Van  Wijngaarden,  which  were  derived  systematically  by  Caflisch,  Miksis,  Papanicolaou  and 
Ting  [23],  We  obtain  numerical  solutions  of  the  nonlinear  system  and  compare  them  with  the  the  linear¬ 
ized  solutions.  We  find  that  even  for  a  small  amplitude  incident  pressure  wave,  it  is  possible  to  have  non¬ 
linear  transmitted  and  reflected  waves.  We  also  consider  the  limit  where  the  bubbly  layer  is  thin  relative 
to  the  incident  wave  length  X. 

By  using  matched  asymptotic  expansions  we  find  that  we  can  replace  the  bubbly  layer  by  an  inter¬ 
face  subjected  to  the  condition  of  continuity  of  pressure  and  a  nonlinear  jump  condition,  i.e. 

P+=P~ 


and 


(u+  -  u  )-n  =  J  KR2R,  dn 

— oo 

in  which  u1  denote  the  velocities  on  the  ±  sides  of  the  thin  bubbly  layer,  n  is  the  unit  normal  vector,  and  n 
and  R  are  the  distance  normal  to  the  layer  and  the  bubble  radius  scaled  by  the  reference  length  h.  Here 
the  parameter  K  is 

„  3(*opc2  h 

K~  P.  X 

where  po  is  the  reference  gas  volume  fraction,  p  and  c  are  the  density  and  the  speed  of  sound  on  the  + 
side,  and  pe  is  the  equilibrium  pressure.  The  equation  for  the  jump  of  the  normal  velocity  contains  inter¬ 
nal  effects  of  the  bubbly  layer.  Solutions  of  this  limiting  case  are  compared  with  the  numerical  results  and 
good  agreement  is  found  even  when  the  ratio  of  the  bubbly  layer  thickness  to  the  wave  length  is  of  order 
one.  These  results  were  presented  at  "A  Workshop  on  Theoretical  Aspects  of  Multiphase  Flow 
Phenomena" ,  at  Cornell  University,  Mathematical  Sciences  Institute,  in  October  1988  and  reported  in 
[24], 


We  also  derived  an  analytical  model  based  on  the  principal  of  minimum  energy  to  determine  the 
steady-state  two-phase  mixture  level  and  axial  gas  volume  fraction  a.  The  gas  is  introduced  uniformly 


across  the  base  of  a  shallow  liquid  pool  and  escapes  through  the  upper  interface  while  the  mass  of  the 
liquid  is  conserved.  The  resulting  Euler  equation  and  the  natural  boundary  condition  at  the  free  interface 
are  reduced  to  algebraic  equations,  i.e. 
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in  which  £  is  the  added  mass  coefficient  It  was  found  that  the  mixture  level  jcj  and  the  gas  volume  frac¬ 
tion  a(x)  are  functions  of  the  Froude  number  Fr  and  the  Weber  number  We,  whereas  the  gas  volume  frac¬ 
tion  oq  just  below  the  free  interface  depends  only  on  the  Weber  number,  i.e. 
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The  gas  volume  fraction  has  an  upper  bound  of  1/2  for  the  limiting  case  of  a  bubbly  homogeneous  flow 
for  which  the  added  mass  coefficient  t,  approaches  infinity.  This  agrees  well  with  experimental  results 
where  bubbly  flows  undergo  a  phase  change  at  an  approximate  void  fraction  of  0.45.  This  investigation 
[25]  will  be  presented  at  the  26th  National  Heat  Transfer  Conference,  August  6-9, 1989,  Philadelphia. 


During  the  period,  September  1987  to  May  1988,  Ting  was  on  sabbatical  leave  from  NYU.  While 
visiting  Stanford  University,  he  and  J.B.  Keller  completed  an  investigation  of  slender  jets  and  thin  sheets 
with  surface  tension  [26].  While  visiting  at  Northwestern  University,  he  and  M.J.  Miksis  completed  a 
study  of  scattering  of  an  incident  wave  from  an  interface  separating  two  fluids  [27], 
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(7)  Helical  Fluid  Flows  (Landman) 

An  important  class  of  fluid  motions  are  those  which  are  constrained  to  be  helical,  i.e.,  in  standard 
cylindrical  polar  coordinates,  the  velocity  of  the  fluid  is  of  the  form  v=v(r,  0+az,r).  By  assuming  this 
symmetry,  the  equations  of  motion  can  be  reduced  to  an  essentially  2-dimensional  polar-coordinate 
representation;  however  the  dynamics  still  allows  the  stretching  of  vortex  lines,  a  3-dimensional 
phenomenon  which  is  of  particular  interest  in  the  study  of  vortex  tube  interactions  and  turbulence. 

A  theoretical  aspect  of  this  work  concerns  the  investigation  of  steady  solutions  of  the  inviscid  Euler 
equations.  In  collaboration  with  Steve  Childress  and  Hank  Strauss  [2],  Landman  has  found  that  the  form 
these  equations  take  allows  for  the  specification  of  two  arbitrary  functions  to  produce  a  steady  helical  vor¬ 
tex.  The  determination  of  these  functions  is  possible  if  one  considers  the  vortex  as  being  formed  in  the 
limit  of  vanishing  viscosity,  known  as  the  Prandtl- Batchelor  limit.  Prandtl  [6]  and  Batchelor  [1]  analyzed 
the  analogous  problem  for  steady  flow  in  2-dimensions  and  with  axisymmetry  (both  are  limits  of  helical 
solutions  as  a  tends  to  zero  and  infinity  respectively).  The  resulting  equations  which  determine  the 
Batchelor  limit  for  the  unknown  functions  are  extremely  complex  for  a  general  helical  geometry  -  they 
are  classified  as  generalized  differential  equations  for  which  there  is  little  theory  available.  The  solution 
of  such  equations  should  produce  helical  vortices  of  relevance  in  the  study  of  propeller  wakes  and  swir¬ 
ling  flows  (such  as  are  produced  in  the  trailing  vortex  of  an  airplane),  and  pipe  flows. 

Landman  has  combined  his  work  on  laminar-turbulent  transition  in  shear  flows  with  the  above 
work,  with  the  idea  that  helical  waves  may  be  the  key  to  understanding  the  nonlinear  development  of 
viscous  flow  in  a  pipe.  Circular  pipe  flow  is  linearly  stable;  however  it  becomes  linearly  unstable  if  the 
pipe  is  rotated  about  its  axis.  He  has  written  a  helical  Navier-Stokes  code,  which  simulates  the  time 
dependent  behavior  of  helical  waves  in  a  pipe.  When  the  pipe  is  rotated,  stable  steady  helical  waves 
which  were  previously  predicted  to  exist  are  found  to  undergo  instabilities.  Quasi-periodic  waves  are 
found,  which  appear  to  undergo  further  bifurcations  to  irregular  (chaotic)  time  dependence  [3].  This 
behavior  seems  to  arise  because  of  the  simultaneous  existence  of  several  interacting  oscillatory  modes 
corresponding  to  those  predicted  by  linear  stability  theory. 

Simulations  performed  when  the  pipe  is  stationary  have  not  revealed  any  finite  amplitude  helical 
waves,  though  these  are  predicted  in  a  high  Reynolds  number  theory  by  Smith  and  Bodonyi  [6].  Further 
simulations  will  be  performed;  however  it  appears  that  for  Reynolds  number  less  than  5000  helical  pipe 
flow  is  absolutely  stable  and  non-axisymmetric  disturbances  decay  on  a  fast  convective  time  scale[4]. 


This  evidence  refutes  any  theory  of  pipe  turbulence  based  on  the  secondary  instability  of  spiraling  waves. 

Nevertheless  it  is  possible  that  transition  in  a  pipe  is  triggered  by  a  small  amount  of  rotation  in  the  inlet 

flow,  implying  that  the  helical  waves  studied  in  the  rotating  pipe,  as  described  above,  are  relevant  in  the 

study  of  stationary  pipe  flow. 
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(8)  Phase  Transitions  and  Their  Stability  for  Mixed  Type  Equations  (Affouf,  Caflisch) 
The  isothermal  fluid  equations  with  viscosity  and  capillarity  terms  are 


v,-ux  =  0 


0) 


ul+p(v)x  =  euxx-Svxxx. 


(2) 


in  which  v  is  specific  volume,  u  is  velocity,  x  is  the  Lagrangian  space  variable,  e  is  the  viscosity 
coefficient  and  6  is  the  capillarity  coefficient.  The  capillarity  term  was  first  proposed  by  Korteweg  [9]  and 
later  analyzed  by  Felderhof  [3]  and  Bongiomo,  Scriven  and  Davis  [2].  This  system  with  a  van  der  Waals 
pressure  law  was  analyzed  by  Slemrod  [15,16]  as  a  simple  model  for  liquid-gas  phase  transitions.  The  van 
der  Waal’s  pressure  term  is 


p(v)  = 


RB 
v  -b 


a 

„2 


(3) 


for  constants  a,  b,  R,  0.  We  use  mainly  the  simplified  form  p(y)  =  v  -  v3.  For  these  choices  of  p,  there 
is  an  interval  [a,  p]  of  values  of  the  specific  volume  v  for  which  p'(v)  >  0.  This  unstable  spinodal  region 
separates  the  two  phases  of  liquid  v  <  a  and  gas  p<v.  For  v  in  the  interval  [a,  p]  the  system  (l)-(3)  is 
elliptic  in  x  and  r,  while  in  the  liquid  or  gas  phases  it  is  hyperbolic. 

The  system  (l)-(3)  is  only  a  crude  model  for  describing  dynamic  phase  transition,  without  a  solid 
physical  basis.  However  we  believe  that  the  mathematical  properties  of  this  system  will  occur  in  other 
more  realistic  phase  transition  models  [4]  or  in  other  conservation  laws  of  mixed  type  [10,14]. 

Traveling  wave  solutions  for  (l)-(3)  representing  phase  transitions  and  shock  waves  were  analyzed 
by  Slemrod  [15,16]  and  Hagan  and  Slemrod  [5],  In  particular  they  showed  that  the  capillarity  term  in  (2) 
is  necessary  to  get  a  full  set  of  phase  transition  waves.  Moreover,  inclusion  of  capillarity  eliminates  some 
phase  jumps  and  thus  reduces  the  multiplicity  of  solutions  for  the  Riemann  problem.  These  traveling 
wave  solutions  were  then  combined  with  rarefaction  waves  by  Shearer  [11-13]  and  Hattori  [6,7]  to  pro¬ 
duce  approximate  solutions  of  the  Riemann  problem.  The  rarefaction  waves  would  be  exact  if  e  and  5 
were  zero.  By  a  Riemann  problem  solution,  we  mean  any  combination  of  traveling  waves  (shocks  or 
phase  transitions)  and  approximate  rarefactions  that  connect  two  given  states  («/,  vt)  at  *=-<»  and  («r,vr) 
at  *=*».  A  surprising  result  of  Shearer  [13]  is  that  the  Riemann  problem  solution  for  given  end  states 
( u/.v, )  and  (ur,vr)  is  not  unique. 
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In  [1]  we  have  performed  numerical  computations  of  the  initial  value  problem  for  (1)  and  (2),  with 
the  pressure  law  (3),  for  rarefaction  waves  and  for  shock  waves  and  phase  transition  waves.  In  particular 
we  numerically  showed  that  all  of  these  waves  are  stable  to  sufficiently  small  initial  perturbations,  if  the 
end  states  of  the  wave  are  not  in  the  spinodal  region.  This  was  surprising  because  of  the  non-uniqueness 
of  Riemann  problem  solutions  and  because  the  phase  transition  waves  pass  smoothly  through  the  spinodal 
region.  Stability  of  these  phase  transitions  can  be  explained  by  their  narrow  width.  On  the  other  hand,  for 
moderate  sized  initial  perturbations,  these  waves  may  be  unstable. 
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VH.2  Numerical  Methods  for  Fluid  Dynamics 


(9)  Convergence  of  the  Point  Vortex  Method  (Goodman,  Hou,  Lowengrub) 

The  vortex  method  for  inviscid,  incompressible  flow  is  an  effective  numerical  method  that  intro¬ 
duces  no  artificial  viscosity.  Convergence  studies  of  the  method  required  smoothing  of  the  vortices  to 
vortex  blobs  with  large  blob  size,  while  point  vortices  or  small  blobs  were  used  in  practice.  Recently 
Goodman,  Hou  and  Lowengrub  [9,10,12]  have  proved  convergence  of  the  point  vortex  method  and  the 
vortex  blob  method  with  smoothing  of  any  size.  This  is  an  unexpected  result  that  is  likely  to  change  our 
thinking  of  the  way  vortex  methods  work. 

The  incompressible,  inviscid  Euler  equations  in  two  dimensions  in  vorticity  form  are 

to,  +  (u-V)m  =  0.  (1) 

where  u  is  the  fluid  velocity  and  to  is  the  vorticity.  The  vortex  method  involves  the  tracking  of  particle 
trajectories.  We  use  the  flow  map,  defined  by 

x  :  R2  x  [0,  T]  R2 


so  that  x(a,  r)  is  the  trajectory  of  the  fluid  particle  which  at  time  t= 0  is  at  the  point  a.  The  trajectory  for 
fixed  a  is  computed  from  the  velocity  field  by 


dxja,  t ) 
dt 


=  u(x(a,  r),  t), 


x(a,  0)  =  a. 


(2) 

(3) 


Incompressibility  implies  that  this  map  is  area  preserving.  The  velocity  is  recovered  from  the  vorticity  by 
the  Biot-Savart  Law 

u  =  K*co=  J  K(x - x') (o(x')dx'  (4) 

R2 


The  Biot-Savart  Kernel  K(x)  has  a  singularity  of  the  form  0  (|x|_d+1 )  at  the  origin,  where  d  is  the  number 
of  space  dimensions.  The  vortex  method  consists  in  finding  approximate  particle  positions  by  discretizing 
(4)  and  solving  the  system  of  ordinary  differential  equations. 
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The  discretization  of  (4),  however,  is  tricky  due  to  the  singularity  of  K.  To  alleviate  this  difficulty, 
Chcrin  in  1973  introduced  the  idea  of  a  cut-off  kernel,  K§,  which  is  non-singular  [7],  Physically,  this 
model,  and  its  subsequent  refinements,  use  blobs  of  vorticity  of  a  certain  shape  and  with  size  8,  a  smooth¬ 
ing  parameter,  to  approximate  the  vorticity  field.  This  is  called  the  vortex  blob  method.  This  method  has 
been  shown  to  be  convergent  by  (e.g.[l,2,3,4,8,l  1,15])  provided  that  the  blobs  are  overlapping.  In  partic¬ 
ular,  if  (4)  is  discretized  by  a  length  scale  h  then  it  is  required  that  5  »  h.  Further,  the  method  is  higher 
order  accurate  depending  on  the  smoothness  of  the  initial  vorticity  and  the  type  of  smoothing  done  to  K. 

In  computational  practice,  one  usually  takes  8  <  h  or  8  =  0  [14],  and  so  the  previous  theory  leaves 
this  computational  practice  in  doubt.  The  point  vortex  method  is  exactly  the  case  when  8  =  0.  It 
corresponds  to  an  approximation  of  the  vorticity  field  by  a  collection  of  point  vortices.  This  method  was 
first  proposed  by  Rosenhead  in  1932  to  simulate  the  motion  of  vortex  sheets  in  two  dimensions  [16].  The 
method  has  the  following  difficulty,  however.  Since  there  is  no  smoothing,  if  two  point  vortices  move 
towards  each  other,  the  velocity  that  each  induces  on  the  other  becomes  unbounded  due  to  the  singularity 
of  K.  The  convergence  of  the  point  vortex  method  and  the  vortex  blob  method  for  vortex  sheets  has  been 
demonstrated  numerically  by  Krasny  [13]  and  analytically,  for  short  times,  by  Caflisch  and  Lowengrub 
[5]  under  the  requirement  of  analyticity.  This  restriction  of  analyticity  seems  to  be  required  because  of  the 
ill-posedness  of  the  vortex  sheet  problem  itself.  In  fact,  it  was  widely  thought  that  the  point  vortex 
method  would  not  converge  in  any  finite  Sobolev  space  even  for  smooth  flows.  But,  recently  for  smooth 
flows,  Goodman,  Hou  and  Lowengrub  have  been  able  to  prove  the  stability  and  convergence  of  the  point 
vortex  method  [9]  and  the  vortex  blob  method  with  small  blobs  (  8  <  h)  [10]  in  two  dimensions.  In  three 
dimensions,  Hou  and  Lowengrub  [12]  proved  convergence  of  a  point  vortex  method.  These  results  are 
valid  for  arbitrarily  long  times. 

The  analysis  is  based  on  the  framework  used  in  the  convergence  study  for  the  vortex  blob  method. 
The  scheme  is  proved  to  be  stable  in  lp  provided  that  the  discretization  error  is  better  than  first  order 
accurate.  Then  by  a  careful  analysis  which  takes  into  account  the  oddness  of  the  Biot-Savart  kernel  K,  the 
method  is  shown  to  be  formally  second  order  accurate.  Thus,  the  method  is  indeed  stable  and  conse¬ 
quently  converges  for  long  times.  It  is  important  to  note  that  this  is  only  a  conditional  stability  result;  los¬ 
ing  accuracy  may  generate  instability. 

One  crucial  observation  that  motivates  the  stability  result  is  the  following:  two  particles  initially 
separated  by  a  distance  h  will  remain  separated  by  a  distance  0(h)  for  later  times  because  of  the 
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smoothness  of  the  flow  map  and  the  inverse  flow  map.  This  suggests  that  if  the  point  vortex  method  is 
better  than  first  order  accurate,  then  any  two  particle  trajectories  are  separated  by  a  distance  0(h)  asymp¬ 
totically.  Therefore,  the  singular  kernel  in  the  point  vortex  method  has  a  natural  cut-off  of  order 
8  =  0(h).  In  this  sense,  the  point  vortex  method  is  equivalent  to  the  blob  method  with  smoothing  size 
8  =  ch.  This  identification  allows  us  to  give  a  particularly  simple  proof  of  convergence  which  relies  on  the 
previous  theory  [12]. 

The  proof  of  convergence  of  the  vortex  method  with  small  blobs  essentially  follows  that  of  the 
point  vortex  method  [10].  This  is  because  the  vortex  method  has  a  natural  cutoff  6  =  0(h)  which  is 
independent  of  blobs  of  size  8  c  h.  Asa  consequence,  smoothing  has  little  effect  if  the  blob  size  is  small 
compared  to  h.  It  is  also  interesting  to  note  that  the  higher  order  accuracy  of  the  blob  method  is  lost  when 
the  blobs  become  0(h)  in  size  or  smaller.  The  method  always  reduces  to  second  order  accuracy  in  this 
case. 


The  convergence  proof  of  the  point  vortex  method  in  three  dimensions  also  roughly  follows  the  two 
dimensional  proof.  There  is  an  additional  difficulty  in  three  dimensions,  though,  which  is  that  the  vorti- 
city  can  be  stretched.  This  problem  can  be  handled  using  the  techniques  of  Beale  and  Majda  [3,4]. 

Finally,  the  theoretical  convergence  results  for  two  dimensional  flows  have  been  confirmed  by 
several  numerical  experiments  [9,10],  The  experiments  indicate  second  order  accuracy,  for  a  fixed  time, 
as  the  numerical  parameters  h  and  8,  in  the  case  of  small  blobs,  are  refined.  The  numerical  experiments 
were  performed  on  a  problem  where  the  exact  solution  to  the  Euler  equations  is  known. 
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(10)  An  Adaptive  Vortex  Method  for  Two-dimensional  Viscous,  Incompressible  Flows 
(Hou,  Papageorgiou,  Papanicolaou,  Zhu) 

J.  Zhu  [ 1  ]  introduced  in  his  Ph.  D.  dissertation  an  adaptive  vortex  method  as  an  alternative  to  the 
circular  vortex,  vortex  sheet  method  for  simulating  numerically  two-dimensional,  vorticity  dominated, 
incompressible  flows  at  h;gh  Reynolds  numbers.  The  basic  vortex  elements  are  ellipses  with  variable 
orientation  and  ellipticity.  They  are  convected  by  the  flow,  viscosity  is  simulated  by  a  random  walk  dis¬ 
placement  and  new  vortex  elements  are  created  at  the  boundaries  to  simulate  vorticity  generation.  In  pre¬ 
vious  implementations  of  vortex  methods,  all  vortex  elements  are  circular  in  the  interior  and  flat  (sheets) 
near  the  boundaries.  By  allowing  the  shape  of  the  vortex  elements  to  change  with  the  flow,  we  have  a 
more  natural  representation  of  the  flow  field,  especially  regarding  its  very  different  structure  near  and 
away  from  boundaries.  T.  Hou,  G.  Papanicolaou  and  D.  Papageorgiou  collaborated  in  this  work  [2],  T. 
Hou  has  done,  in  addition,  a  theoretical  analysis  of  a  version  of  the  adaptive  vortex  method  [3]. 

Consider  for  example  a  boundary  layer  flow  described  by  the  Prandtl  equations.  They  are  derived 
after  stretching  the  coordinate  perpendicular  to  the  boundary  and  neglecting  terms  of  order  1  /^Re  or 

smaller.  Suppose  the  local  boundary  is  approximated  by  '■>  =  0  and  the  vorticity  0)  =  ~  -  r—-  by  — r^-. 

dx  dy  dy 

The  vortex  method  implementation  in  this  approximation  is  the  vortex  sheet  method,  which  works  well  in 
the  regime  where  the  Prandtl  equations  hold.  But  the  boundary  layer  is  only  part  of  the  flow  and  the  rest 
has  to  be  obtained  using  the  Navier-Stokes  equations,  which  can  be  done  by  the  circular  vortex  biob 
method,  rhis  requires  matching  flows  in  the  two  regimes  to  generate  the  global  solution.  Since  the  nature 
of  a  vortex  sheet  is  quite  different  from  that  of  a  vortex  blob,  the  transition  in  the  flow  from  the  boundary 
layer  solution  to  the  solution  in  the  interior  is  not  smooth.  With  the  elliptic  vortex  element  implementa¬ 
tion  we  eliminate  this  sudden  transition  by  employing  a  single  solution  which  is  good  in  both  regions. 

Based  on  »he  Prandtl  boundary  layer  approximation,  the  idea  of  the  elliptic  vortex  method  is  this. 

If  we  use  the  same  stretching  as  that  in  the  derivation  of  the  Prandtl  equations  and  keep  the  term  —  in 

dx 

the  vorticity  formulation,  the  vortex  element  is  going  to  be  an  elliptic  blob  with  an  order  ellipticity. 
This  is  the  fixed  elliptic  vortex  blob  method  associated  with  the  transformation  that  is  the  stretching  in  the 
direction  perpendicular  to  the  boundary.  We  extend  this  idea  to  general  continuously  differentiable 
transformations  over  the  flow  region  to  realize  a  smooth  transition  between  interior  and  boundary  layer. 
Given  a  time-dependent  transformation  <J>(x,y,r)  =  (<j>,\y),  the  flow  region  is  mapped  to  a  region  in  a  new 
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coordinate  system  0*,y*).  We  write  the  fluid  equations  in  this  transformed  space  and  define  trajectories 
of  fluid  particles,  velocity  and  vorticity.  We  formulate  a  vortex  method  in  this  transformed  system  assum¬ 
ing  that  all  the  blobs  in  the  new  system  are  circular  and  that  the  vorticity  is  constant  inside  each  blob.  We 
then  discretize  the  kernel  which  corresponds  to  the  Biot-Savart  kernel  in  the  original  system.  Suppose  the 
blob  size  is  so  small  that  linearization  of  O  is  allowed  in  the  neighborhood  of  each  blob.  Then  each  circu¬ 
lar  blob  in  the  new  system  will  represent  an  elliptic  blob  in  the  original  system.  If  we  have  a  collection  of 
circular  blobs  in  the  transformed  space,  the  transformation  d>  will  furnish  us  with  a  collection  of  elliptic 
vortex  blobs  in  the  original  space,  each  determined  by  the  linear  transformation 

x*  =  x*  +0,(x(,r)(x-x1) 

associated  with  the  center  of  the  vortex  blob,  where  x,  and  x*  are  the  positions  of  rth  blob  in  original  and 
transformed  space  respectively,  and  <t>x  is  the  Jacobian  matrix.  It  remains  therefore  to  solve  for  a  collec¬ 
tion  of  elliptic  vortex  blobs  with  constant  vorticity  inside  each  blob.  We  have  the  analytic  solution  for  the 
flow  induced  by  each  elliptic  vortex  blob  from  potential  theory  and  since  the  equation  relating  vorticity  to 
velocity  is  linear,  the  velocity  field  associated  with  the  collection  of  vortex  blobs  is  known  everywhere. 
The  updating  of  vorticity  is  done  in  the  same  way  as  in  the  standard  method.  We  move  each  vortex  blob 
according  to  the  local  velocity  and  a  random  walk  step  is  added  to  simulate  the  diffusion  part. 

This  formulation  bridges  the  gap  between  the  vortex  sheet  method  and  the  circular  vortex  blob 
method  in  a  smooth  way,  with  the  full  system  of  Navier-Stokes  equations  used  in  all  steps.  We  are 
required  to  provide  a  transformation  <t>  adapted  to  the  flow  in  some  way.  Since  this  transformation  is  an 
auxiliary  tool  designed  to  insure  that  the  correct  equations  are  solved  while  implementing  the  vortex 
method  with  as  much  flexibility  as  possible,  the  choice  of  <t>  is  virtually  unrestricted.  We  have  proposed 
two  criteria  for  based  on  different  numerical  considerations.  The  one  we  use  in  our  calculations  is 
obtained  from  the  real  flow  by  averaging  and  truncation  and  satisfies  the  same  boundary  condition  as  the 
real  flow.  The  other,  which  we  formulated  but  did  not  use,  determines  the  transformation  by  a  global 
minimization  of  the  rate  of  shear  in  each  time  step. 

We  performed  several  numerical  experiments  to  test  the  adaptive  vortex  method.  First  we  applied  it 
to  the  flow  past  a  flat  plate.  Our  results  compare  well  with  the  steady  Blasius  solution  after  a  sufficiently 
long  time.  We  also  applied  this  method  to  the  unsteady  flow  in  a  channel  with  a  back  step.  We  compared 
our  results  both  with  experimental  results  and  with  numerical  results  using  the  vortex-blob-sheet  method 
as  performed  by  Setluan  and  Ghoniem  [4],  Good  qualitative  agreement  is  obtained  in  both  comparisons. 


The  adaptive  vortex  method  has  several  advantages.  It  does  not  use  the  Prandtl  equations,  so  that 
we  can  apply  it  to  problems  with  curved  boundaries  without  having  to  prescribe  a  priori  length  scales 
associated  with  the  boundary  layers,  and  sharp  transitions  in  the  flow  are  caught  in  detail.  It  also  avoids 
using  different  sets  of  equations  in  different  regions.  We  have  instead  a  continuously  varying  vortex- 
element  representation  suitable  to  specific  regions.  Since  we  make  no  simplification  of  the  Navier-Stokes 
equations  we  do  not  have  to  worry  about  applicability  of  the  method  to  difficult  problems  where 
backflow  and  separation  appear.  This  is  evident  in  the  backstep  flow  problem,  since  the  vortex  sheets  are 
not  very  good  at  representing  the  flow  in  the  step  comer  where  the  local  Reynolds  number  can  be  small 
and  the  Prandtl  approximation  therefore  inappropriate. 
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(11)  A  Particle  Method  for  the  Numerical  Computation  of  Dilute  Fluid-Particle  Flows 
(Papanicolaou,  Zhu) 

The  continuum  modeling  of  multiphase  flows  has  received  a  lot  of  attention  by  theoreticians,  but  it 
is  generally  agreed  that  only  in  very  special  situations  are  the  continuum  equations  well-founded  and  can 
be  used  with  confidence.  One  such  case  is  fluid-particle,  two-phase  flow  with  the  particle  phase  a  dilute 
or  low  volume-fraction  component.  Dilute  suspensions  of  particles  in  fluids  or  slightly  dusty  gases  are 
good  examples.  A  very  good  survey  of  theoretical  developments  in  dusty  gases  is  Marble’s  paper  [1  ]. 

Papanicolaou  and  Zhu  [2,3]  have  formulated  a  continuum  model  for  two-dimensional  fluid-particle 
flow  in  which  the  fluid  is  incompressible  and  nearly  inviscid  (high  Reynolds  number  relative  to  macros¬ 
copic  length  scales)  while  the  particle  phase  is  dilute  and  the  fluid-particle  interaction  is  calculated  by  the 
Oseen  formula.  We  assume  therefore  that  the  Reynolds  number  relative  to  the  size  of  the  particles  is 
small.  There  are  no  boundaries  in  the  flow  region. 

To  compute  such  flows  numerically  we  use  a  vortex-blob  method  to  represent  the  fluid  phase  and  a 
particle  method  to  represent  the  particle  phase.  The  physical  particles  and  the  “numerical”  particles  are 
very  distinct  entities,  as  is  always  the  case  with  particle  methods  in  numerical  computation.  Modeling  of 
the  interaction  between  the  two  phases  at  the  numerical  level  is  the  main  difficulty  that  must  be  faced  in 
the  problem.  We  have  done  extensive  numerical  experiments  to  test  the  performance  of  various  interac- 
i  tion  models.  We  have  found  that  the  use  of  Voronoi  diagrams  both  in  the  fluid  and  particle  phase  is  very 
effective.  The  necessary  interpolations  are  done  by  an  algorithm  [4]  that  is  suitably  local.  Extensive 
numerical  testing  has  shown  that  our  vortex-particle  method  produces  good  results  for  flows  in  the  regime 
I  described  above. 

Previous  computations  of  fluid-particle  flows  at  high  Reynolds  numbers  using  the  vortex  method 
[5,6]  treated  the  particle  phase  as  tracer  particles  advected  by  the  fluid  but  not  affecting  it  by  their  motion. 

! 
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(12)  Lagrangian  Finite  Element  Method  (Chacon,  Hou) 

Tomas  Chacon  and  Thomas  Hou  [1]  have  proposed  a  new  Lagrangian  finite  element  method  to 
approximate  the  2-D  incompressible  Euler  equations,  in  this  method,  they  transport  the  initial  triangula¬ 
tion  by  streamlines  of  the  flow  and  approximate  the  vorticity  by  piecewise  linear  finite  elements.  As  a 
consequence  the  vorticity  is  exactly  computed  on  the  mesh  points  of  the  tri angulation.  Their  method  is 
conceptually  similar  to  the  vortex  method,  except  that  they  do  not  need  to  introduce  smoothing  for  the 
singular  kernel.  This  is  possible  because  the  product  of  the  Biot-Savart  kernel  and  a  polynomial  can  be 
integrated  analytically  over  any  triangle.  This  enables  them  to  obtain  a  stable  and  consistent  discretization 
for  the  singular  velocity  integral  without  smoothing.  The  kind  of  stability  result  they  obtain  here  is  a  con¬ 
ditional  one  in  that  the  method  is  stable  if  the  consistency  error  is  better  than  first  order  accurate.  In  order 
to  ensure  stability,  they  need  to  use  linear  or  higher  order  finite  element  approximation  for  the  vorticity 
field  and  a  better  than  first  order  time  integration  formula.  Then  the  Lagrangian  triangulation  will  not 
become  degenerate  provided  h  <ho(T)  is  small.  Consequently  the  method  is  stable  and  convergent. 
Numerical  experiments  have  shown  that  the  method  converges  with  second  order  accuracy  for  long  times. 

The  analytical  techniques  used  here  are  of  independent  interest,  since  they  do  not  need  to  differen¬ 
tiate  the  singular  kernel  in  proving  the  stability  as  is  the  case  of  the  vortex  method.  Only  the  smoothness 
of  the  solution  and  the  L\oc  property  of  the  kernel  are  used  in  the  stability  argument.  This  enables  them  to 
prove  convergence  of  the  method  for  a  more  singular  kernel  than  that  of  the  2-D  Euler  equations  [3). 
Another  application  of  this  method  is  in  calculating  the  motion  of  vortex  patches.  The  fact  that  the  kernel 
is  integrated  exactly  over  each  triangle  makes  the  finite  element  method  very  natural  for  approximating 
vortex  patchs.  Numerical  computations  have  shown  that  the  method  provides  accurate  approximations  for 
two  interacting  elliptical  patches  up  to  the  time  when  the  patches  are  nearly  singular.  A  convergence 
analysis  for  vortex  patch  problems  is  in  progress  [2]. 
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VTL3  Random  Media  and  Composite  Materials 

(13)  Direct  and  Inverse  Problems  for  Pulse  Reflection  of  Waves  from  Random  Media 
(Asch,  Kohler,  Papanicolaou,  Postel,  Sheng,  White) 

Consider  a  pulsed  acoustic,  electromagnetic  or  elastic  wave  impinging  on  a  material  that  is  inho¬ 
mogeneous,  the  earth’s  crust  for  example.  The  reflected  signal  is  measured  at  the  surface  and  it  looks  very 
noisy,  like  a  seismogram.  What  information  about  the  medium  can  we  extract  from  this  signal?  In 
reflection  seismology  this  is  a  classic  question  that  has  received  a  great  deal  of  attention  over  the  last  forty 
years.  It  has  also  received  a  lot  of  attention  as  a  difficult  problem  in  mathematics:  the  inverse  scattering 
problem.  George  Papanicolaou  in  collaboration  with  R.  Burridge  (now  at  Schlumberger-Doll  Research), 
P.  Sheng  and  B.  White  (at  Exxon  Research  and  Engineering  Co.)  and  more  recently  W.  Kohler  (Virginia 
Polytechnic  Inst.),  M.  Postel  (Courant  Institute)  and  M.  Asch  (a  graduate  student  at  the  Courant  Institute) 
began  a  study  of  this  question  as  a  problem  in  stochastic  processes.  We  model  the  reflecting  material  as  a 
random  medium,  characterize  all  possible  reflected  signals  as  a  particular  class  of  nonstationary  stochastic 
processes  and  then  pose  the  recovery  of  material  properties  as  a  statistical  estimation  problem  for  this 
class  of  processes.  We  test  the  theory  by  generating  numerically  reflected  signals  and  then  using  our  esti¬ 
mation  algorithm  on  them.  References  [1-10]  contain  our  work.  Reference  [11]  is  a  recent  textbook  on 
seismology  where  a  variety  of  inverse  problems  are  posed  and  references  [12-13]  are  the  papers  that  gave 
us  the  idea  to  start  this  research.  They  contain  numerical  simulations  for  pulse  reflection  from  random 
media  but  they  do  not  address  the  inversion  questions. 

How  can  we  develop  a  theory  for  reflected  signals  that  is  realistic  but  simple  enough  to  be  useful? 
First  we  must  choose  carefully  a  model  for  the  reflecting  medium  and  then  analyze  its  scattering  proper¬ 
ties  with  methods  that  are  well  suited  to  the  model.  We  have  focused  attention  on  layered  acoustic  half 
spaces  so  far.  Locally  layered  elastic  half  spaces  can  also  be  studied  by  our  methods  and  they  are  rather 
realistic  models  for  the  earth.  The  next  issue  is  modeling  of  the  inhomogeneities,  which  means  that  we 
must  decide  what  form  the  density  and  bulk  modulus  of  the  material  are  to  have.  We  have  already 
assumed  that  they  are  functions  of  one  variable  only:  a  layered  medium.  Motivated  by  experience  with 
other  problems  in  materials  science,  we  decided  to  model  the  medium  properties  by  random  functions 
that  vary  on  two  length  scales.  Their  mean  values  vary  on  a  macroscopic  scale,  for  example  kilometers  in 
the  earth  model.  Their  fluctuations  are  rapidly  varying  but  not  necessarily  small.  In  the  earth  model  the 
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speed  of  sound  varies  on  the  scale  of  two  or  three  meters  by  as  much  as  thirty  percent,  as  we  know  from 
well  log  data.  Of  course  there  are  variations  on  many  other  scales  in  the  earth  but  as  long  as  there  are  two 
well  separated  scales  that  predominate,  our  methods  of  analysis  apply.  The  final  step  in  the  modeling  is 
the  choice  of  the  probing  pulse  that  will  generate  the  reflected  signal.  We  choose  the  pulse  width  to  be  on 
an  intermediate  scale  between  the  macro  and  microscales.  In  the  earth  model  the  pulse  width  may  be  150 
meters  (a  20-40  Hz  pulse)  which  is  narrow  compared  to  the  10  kilometer  depth  we  want  to  probe  but 
broad  compared  to  the  3  meter  width  of  a  typical  layer. 

Now  the  problem  is  formulated  as  one  in  the  asymptotic  analysis  of  stochastic  equations  on  which 
we  have  worked  extensively  in  the  past  and  have  developed  mathematical  techniques.  The  main  outcome 
of  the  asymptotic  theory  is  a  precise  characterization  of  the  reflected  signals  as  a  relatively  simple  class  of 
nonstationary  Gaussian  processes  that  depend  primarily  on  the  macroscopic  properties  of  the  medium  and 
are  practically  independent  of  the  small  scale  variations. 

Once  the  scope  of  the  theory  is  fully  understood  we  can  go  to  the  inverse  question  which  is  an 
unusual  estimation  problem  for  a  Gaussian  process.  The  main  advantage  of  our  theory  is  that  we  recover 
directly  the  large  scale  properties  of  the  medium  and  never  get  tangled  with  reconstructing  all  its  details. 
This  is  the  main  drawback  of  other  methods,  overcome  mainly  by  ad  hoc  "processing"  of  the  reflected 
signals. 

To  implement  the  theory  we  must  first  have  an  efficient  way  to  generate  numerically  reflected  sig¬ 
nals  for  our  model.  This  has  been  done  for  the  earth  model  with  the  parameters  given  above  when  a  plane 
wave  pulse  is  incident.  A  finite  difference  method  is  used.  To  resolve  the  small  scales  (3  meters)  we  use  a 
spatial  mesh  that  is  about  30  cm.  Then  we  solve  for  long  enough  times  that  a  signal  can  go  down  8-10 
kilometers  and  return.  This  takes  about  150,000  time  iterations  when  the  speed  is  3  km/sec  on  average.  A 
well-vectorized  code  on  the  cyber  205  produces  a  reflected  signal  with  these  parameters  in  about  25  cpu 
minutes. 

Once  we  have  good  reflected  signals  (synthetic  seismograms)  we  can  try  our  inversion  algorithm. 
This  was  done  in  detail  in  [4]  and  the  results  are  quite  good  if  enough  redundancy  is  used.  That  is,  if  we 
use  enough  independent  samples  of  reflected  signals  in  our  estimation  algorithm.  The  quantity  estimated 
is  the  smooth  (averaged)  part  of  the  local  sound  speed.  For  a  good  estimation  we  need  about  50  samples. 
But  in  reality  we  have  only  one  reflected  signal.  How  do  we  generate  other  realizations?  In  reflection 
seismology  the  exciting  pulse  is  not  a  plane  wave;  it  is  also  localized  in  space  because  it  comes  from  a 

I 


-46- 


pulsed  point  source.  We  then  measure  reflected  signals  at  different  points  on  the  surface  and  treat  them  as 
separate  realizations.  This  must  of  course  be  done  correctly  and  the  theory  to  do  it  is  in  [5].  The  main 
point  of  our  computational  research  in  the  immediate  future  is  the  numerical  simulation  of  such  reflected 
signals. 

Generating  reflected  signals  from  a  pulsed  point  source  over  a  randomly  layered  medium  is  a  very 
big  computational  task.  For  a  plane  wave  we  used  a  spatial  resolution  of  30  cm  in  a  finite  difference 
scheme.  In  a  two  or  three  dimensional  setting  such  resolution  is  out  of  the  question  when  we  want  propa¬ 
gation  over  kilometers.  In  dimensions  higher  than  one  we  may  use  a  layer-specific  method  for  computing 
the  wave  without  having  to  resolve  it  inside  the  layers.  This  requires  taking  Fourier  transforms  in  space 
horizontally  and  in  time  with  enough  resolution  to  capture  the  space-time  fluctuations  in  the  signal.  If  we 
use  2000  layers  (6  km  maximum  depth)  it  becomes  clear  that  this  is  again  a  very  big  problem,  even  for 
generating  one  realization.  We  have  at  present  a  well-vectorized  code  for  this  problem  (layer  algorithm 
plus  Fourier  and  Fourier-Hankel  transforms)  running  quite  efficiently  on  an  ETA  10.  It  is  producing  very 
accurate  results  but  it  still  takes  about  30  cpu  minutes  on  an  ETA  10  to  generate  one-second  signals  at 
about  20  surface  locations.  We  are  also  developing  on  an  IRIS  4D/50  GT  the  graphics  necessary  to  visu¬ 
alize  the  results.  We  would  like  to  be  able  to  generate  3-4  second  reflected  signals  at  50  to  100  points  on 
the  surface  in  different  spatial  configurations  so  that  we  can  test  our  theory  and,  more  important,  test  vari¬ 
ous  inversion  schemes. 

W.  Kohler,  G.  Papanicolaou  and  B.  White  are  completing  a  monograph  on  the  theory  that  has  been 
developed,  which  will  be  published  by  North  Holland  in  1989. 
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(14)  Reflection  and  Transmission  of  Gravity  Waves  in  a  Channel  with  Irregular  Bottom  Topogra¬ 
phy  (Nachbin,  Papanicolaou) 

The  propagation  of  long  gravity  waves  in  a  channel  can  be  analyzed  using  the  shallow  water  equa¬ 
tions.  For  weakly  nonlinear  waves  traveling  in  one  direction  this  is  the  Korteweg-de Vries  equation.  For 
fully  nonlinear  waves  traveling  in  both  directions  one  can  use  conservation  laws  analogous  to  those  in  gas 
dynamics,  having  a  barotropic  equation  of  state  with  y=2.  For  small  amplitude,  linear  waves,  the  shallow 
water  equations  reduce  to  the  acoustic  equations  [1]. 

The  shallow  water  equations  are  only  valid,  however,  when  the  variations  in  the  bottom  topography 
are  slow  compared  to  the  average  depth.  In  [2]  we  showed  that  if  the  bottom  topography  varies  periodi¬ 
cally  with  the  period  comparable  to  the  depth  and  if  weakly  nonlinear,  long  waves  propagate  in  one  direc¬ 
tion,  then  again  the  Korteweg-de  Vries  equation  is  a  valid  approximation.  The  coefficients  in  it  have  to  be 
determined  by  solving  cell  problems  as  is  usual  in  homogenization  [3).  The  analysis  of  [2]  does  not  gen¬ 
eralize  to  bottom  topographies  that  vary  randomly.  We  expect  substantial  reflection  in  the  random  case, 
which  means  that  the  assumption  of  waves  moving  in  one  direction,  underlying  the  Korteweg-de  Vries 
approximation,  is  no  longer  valid.  Furthermore,  even  in  the  linear  case  there  is  no  clear  understanding  of 
how  to  treat  reflection  and  transmission  in  the  random  case  [4-7]. 

A.  Nachbin  and  G.  Papanicolaou  [8,9]  have  begun  a  systematic  analysis  of  the  reflection- 
transmission  problem  with  random  bottom  topography.  The  analysis  is  in  two  parts  and  is  restricted  to 
linear  problems  at  present  The  first  part  is  the  development  of  an  efficient  code  for  calculating  numeri¬ 
cally  the  reflection  and  transmission  of  wave  pulses  by  random  bottom  inhomogeneities.  We  use  a  boun¬ 
dary  element  method  as  has  been  done  previously  by  several  investigators  [10-16].  In  previous  imple¬ 
mentations  the  boundary  element  method  has  rendered  dispersive  effects  better  than  finite  difference 
methods  of  comparable  complexity.  There  is  a  simple  theoretical  reason  for  this  that  we  have  clarified. 
Apart  from  that  our  implementation  of  the  boundary  integral  method  differs  from  previous  ones  only  in 
that  our  code  is  fully  vectorized  and  runs  at  the  equivalent  of  about  400  megaflops  on  the  ETA  10  at  the 
John  von  Neumann  Center  at  Princeton.  This  allows  us  to  study  propagation  of  pulses  for  a  long  enough 
time  that  multiple  scattering  from  the  bottom  irregularities  is  important.  It  is  also  very  important  for  the 
next  step  in  our  study,  the  study  of  nonlinear  effects. 

The  second  part  of  our  study  is  theoretical.  We  use  the  asymptotic  theory  of  stochastic  equations  to 
analyze  reflection  and  transmission.  We  have  shown  that  although  the  acoustic  approximation  is  not  valid 
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for  long  waves  over  an  irregular  bottom,  the  statistical  description  of  the  reflection-transmission  process 
has  an  equivalent  acoustic  interpretation,  with  the  dispersive  effects  accounted  correctly.  An  important 
feature  of  our  analysis,  which  parallels  that  of  [17],  is  that  the  bottom  inhomogeneities  need  not  be  small. 
It  also  allows  us  to  address  inverse  problems:  can  the  reflected  surface  waves  tell  us  what  the  large-scale 
features  of  the  bottom  topography  are?  This  question  has  a  rather  complete  answer.  Naturally  the  addi¬ 
tional  length  scale  present  in  this  problem,  the  average  depth,  not  present  in  the  acoustic  case,  compli¬ 
cated  the  implementation  of  the  theory  considerably. 
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(15)  Bounds  on  the  Effective  Properties  of  Composite  Materials  (Avellaneda,  Kohn,  Milton) 

Despite  a  surge  of  interest  in  composite  materials  it  is  typically  not  known  how  to  calculate  the 
range  of  possible  effective  conductivity  tensors  as  the  microstructure  is  varied  over  all  configurations. 
This  question  is  relevant  to  optimization.  For  isotropic  poly-crystalline  materials  constructed  from  a 
given  single  uniaxial  crystal,  Schulgasser  [1]  made  a  conjecture  in  1982  that  the  most  resistive 
configuration  would  be  a  sphere  assemblage  in  which  the  crystal  axis  is  directed  radially  outwards  in  each 
sphere.  Avellaneda,  Cherkaev,  Lurie  and  Milton  [2]  have  now  proved  this  conjecture,  using  a  technique 
which  involves  null-Lagrangians.  The  technique  is  evidently  powerful  since  it  also  enabled  us  to  prove  a 
sharp  phase  interchange  inequality  which  had  been  conjectured  in  1982  [3].  This  inequality 

o*o*  /Oj  o2  +  (o*  +  o*  )/(oi  +  cr2)  >  2 

correlates  the  effective  conductivity  o*  of  an  isotropic  composite,  made  from  two  isotropic  phases  of  con¬ 
ductivities  Oi  and  o2  with  the  effective  conductivity  o*  having  the  same  microstructure  but  with  the 
phases  interchanged.  This  inequality  has  several  important  practical  implications  [4]. 

Another  technique  to  bound  the  range  of  possible  effective  tensors  is  to  use  Hashin-Shtrikman  vari¬ 
ational  principles.  Milton  and  Kohn  [5]  have  streamlined  this  method  and  used  it  to  obtain  bounds  on  the 
elasticity  tensor  of  two-component  materials.  The  resulting  bounds  are  natural  generalizations  to  aniso¬ 
tropic  composites,  of  the  widely  known  and  practically  important  Hashin-Shtrikman  bounds.  The  bounds 
are  optimal  and  are  attained  with  the  class  of  sequentially  layered  laminate  materials.  Such  materials  are 
constructed  by  layering  together  the  components  on  widely  separated  length  scales. 

Milton  has  also  investigated  [6]  the  problem  of  conduction  in  two-dimensional,  two  component 
materials  in  which  the  conductivity  tensor  is  asymmetric,  as  it  will  be  in  the  presence  of  a  magnetic  field. 
By  the  use  of  a  suitable  transformation  he  finds  this  problem  is  isomorphic  to  one  in  which  both  com¬ 
ponents  have  symmetric  conductivity  tensors,  i.e.  in  which  the  magnetic  field  is  absent.  This  surprising 
result  enables  him  to  completely  characterize  the  set  of  all  possible  asymmetric  effective  tensors  of  the 
composite  as  the  microstructure  is  varied. 

Milton  plans  to  continue  characterizing  the  set  of  possible  effective  conductivity  tensors  associated 
with  composites  of  various  given  materials.  The  problem  bears  some  similarity  with  the  study  of  phase 
transitions  in  Statistical  Physics.  We  need  to  identify  microstructures,  i.e.  phases,  that  minimize  or  max¬ 
imize  some  linear  combination  of  energies.  These  energies  may  for  example  represent  the  electrical 
energy  stored  in  the  composite  when  the  field  is  applied  in  several  different  test  directions.  In  a  typical 
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problem  "microstructure"  transitions  need  to  be  identified,  and  some  physical  explanation  for  the  transi¬ 
tion  needs  to  be  given.  This  study  is  still  in  its  infancy.  No  one  yet  knows  what  are  the  possible  types  of 
transitions  that  can  occur. 

The  available  methods  for  bounding  effective  constants,  all  seem  to  have  their  limitations.  In  par¬ 
ticular  the  methods  fail  to  give  a  complete  characterization  of  the  set  of  possible  effective  tensors  when  an 
isotropic  material  is  mixed  with  an  anisotropic  material  in  prescribed  proportions  even  in  two- 
dimensions.  The  problem  is  generic  and  stems  from  a  lack  of  knowledge  of  the  extent  to  which  electric 
fields  can  be  concentrated  in  small  volumes  in  the  composite.  A  related  question  is  to  investigate  the  pro¬ 
bability  distribution  of  field  intensities  in  a  composite.  While  this  is  an  intractable  problem  for  most 
materials,  one  can  hope  to  obtain  results  for  a  special  class  of  hierarchical  models.  These  models  are 
interesting  because  they  exactly  realize  [7,8]  the  well-known  effective  medium  approximation  for  the 
conductivity. 

It  may  also  be  possible  to  further  explore  the  connection  between  various  analytic  representations 
for  the  effective  conductivity  as  a  function  of  the  component  conductivities  in  a  multi -component  compo¬ 
site.  From  the  analytic  properties  of  this  function  Golden  and  Papanicolaou  [9]  have  deduced  an  integral 
representation  for  the  effective  conductivity  function.  It  now  seems  their  representation  can  be  modified 
to  incorporate  the  natural  symmetries  amongst  the  roles  of  the  components.  This  should  lead  to  a  very 
elegant  representat  on  formula.  On  the  other  hand  there  is  a  continued  fraction  representation  for  the 
effective  conductivity  [10]  that  was  obtained  by  recursion  methods  in  the  Hilbert  space  in  which  the  con¬ 
ductivity  problem  is  set.  It  remains  to  connect  these  representations  and  to  find  the  analytical  significance 
of  the  parameters  that  enter  the  continued  fraction  expansion.  This  problem,  if  solved,  would  likely  have 
a  profound  influence  in  the  study  of  functions  of  several  complex  variables. 
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(16)  Optimality  of  Bounds  for  the  Polycrystalline  Problem  (Milton,  Nesi) 

Composite  materials  have  received  a  lot  of  mathematical  attention  in  the  last  few  years.  This  is  in 
part  due  to  their  relation  to  the  calculus  of  variations,  optimal  control  theory  and  homogenization  of  fam¬ 
ilies  of  elliptic  p.d.e.’s  with  rapidly  oscillating  coefficients. 

A  particular  goal  of  much  recent  work  is  the  so  called  "G-closure  problem",  which  seeks  to  deter¬ 
mine  the  set  of  all  possible  composites  attainable  from  a  given  class  of  component  materials.  The  thesis 
of  Nesi  [12]  is  concerned  with  this  problem  in  the  context  of  polycrystalline  composites,  for  electric  con¬ 
ductivity  in  three  space  dimensions. 

Polycrystals  are  composite  materials  consisting  of  individual  grains  (or  crystallites)  having  the 
same  basic  conductivity  tensor  but  varying  in  orientation  from  one  to  the  next  The  effective  conductivity 
of  such  a  material  has  been  investigated  by  many  authors  including  [1-6]. 

From  the  mathematical  viewpoint,  one  can  write  the  field  equations  for  the  electric  field‘d and  the 
current  field  y^in  a  polycrystalline  composite  in  the  form 

y(x*)  =  a(x)-?(?),  l?Qt)  =  -grad  u (**),  div~f<ft)  =  0 

where  o(3?)  is  the  local  conductivity  tensor  and  u(jt)  is  the  electrostatic  potential.  By  the  definition  of  a 
polycrystalline  composite,  has  the  form 

o(?)  =  RT(5t)oR  Ot), 

where  a  =  diag{a o2,  Oj}  is  nonnegative  definite  and  represents  the  conductivity  tensor  of  the  basic 
crystallite  and  R(^t)  is  a  field  of  orthogonal  matrices  in  R 3  reflecting  the  variation  of  crystal  orientation  in 
the  phase  geometry.  One  way  to  define  the  effective  conductivity  tensor  o*  is 

~tfTa*l?=  lim  inf  — f(  V  uQt))T  o(-)  V  u(r)  d? 
e->o  u(?)  \D  \  JD  e 

on  dD. 

Significant  progress  in  understanding  the  effective  moduli  of  polycrystalline  composites  has 
recently  been  obtained  in  [7],  That  work  derives  improved  bounds  on  the  effective  conductivity  using  a 
method  based  on  null  Lagrangians,  c.f.  [8]. 

The  goal  of  Nesi’s  work  is  to  explore  the  attainability  of  a  certain  "lower  bound"  proved  in  [7]  for 
the  effective  conductivity  of  a  polycrystalline  composite.  His  main  result  is  that  this  lower  bound  is 
attained  by  a  large  class  of  materials.  When  specialized  to  isotropic  effective  conductivities,  his  result  is 
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as  follows:  given  a  basic  crystal  with  conductivity  tensor  o=diag{<3\,  02,  03J,  the  lowest  possible  isotro¬ 
pic  effective  conductivity  of  a  polycrystalline  composite  made  from  this  material  is 

a*=as  (1) 

where  os  is  the  unique  positive  root  of  the  cubic  equation 

a?+a?(Oi  -KJ2+CJ3MO!  o2o3  =0.  (2) 

His  other  results  are,  in  a  sense,  extensions  of  (2)  to  the  case  of  anisotropic  effective  tensors.  He  proves 
that  the  bound 

a^+0?(0,*-KJ2*+C3*)-4CTi*02*03*  <0.  (3) 

is  saturated  (i.e.  equality  holds)  for  a  large  class  of  "sequentially  laminated"  materials,  corresponding  to  a 
variety  of  different  anisotropic  effective  conductivity  tensors.  Before  his  work  the  attainability  of  (1),  was 
known  only  in  the  very  special  case  when  the  basic  crystal  has  uniaxial  symmetry  (i.e.  oj  =o2  or 
02  =  03).  Nothing  was  known  in  the  general  case  Oi  <  o2  <  03 . 

To  obtain  such  a  result,  one  must  find  a  microgeometry  for  which  the  associated  effective  conduc¬ 
tivity  achieves  equality  in  the  bound  (3).  He  has  found  several  different  constructions  which  do  this.  All 
are  based  on  sequential  lamination.  This  is  a  very  flexible  method  for  constructing  microgeometries  with 
explicitly  computable  effective  moduli;  it  has  been  successfully  applied  to  many  G-closure  problems,  see 
e.g.[ll].  The  problem  with  the  method  of  sequential  lamination  is  that  there  are  too  many  parameters  to 
be  chosen.  In  a  certain  sense,  Nesi’s  contribution  is  to  show  what  are  the  extremal  ways  of  selecting  these 
parameters  in  this  problem. 

The  first  step,  initiated  by  G.W.  Milton,  was  a  set  of  numerical  experiments:  they  looked  for  an  iso¬ 
tropic  material  with  the  lowest  possible  conductivity  by  laminating  the  basic  crystal  with  some  rotation 
of  itself.  By  examining  the  extremal  microgeometries  generated  numerically,  they  deduced  which  stra¬ 
tegy  of  lamination  seemed  to  be  the  most  efficient.  Then  Nesi  proved  analytically  the  attainability  of  the 
lower  bound  (3)  by  introducing  an  increasing  number  of  length  scales  (in  the  context  of  reiterated  homo¬ 
genization,  see  [9]  and  [10])  and  making  use  of  a  classical  fixed  point  theorem.  The  constructions  Nesi 
introduced  give  rise  to  microgeometries  which  were  never  considered  before  in  the  literature.  This  type  of 
microgeometry  is  called  an  infinite  rank  composite  material,  since  the  homogenization  has  to  be  reiterated 
infinitely  many  times  to  get  the  extremal  composite  which  attains  the  bound. 
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Nesi  has  also  studied  schemes  which  correspond  to  finite  rank  lamination.  Such  schemes  depend  on 
a  smaller  number  of  parameters  which,  as  they  vary,  lead  to  trajectories  in  the  space  of  effective  tensors, 
lying  on  the  boundary  of  the  G-closure.  He  has  obtained  a  complete  characterization  of  those  lamina¬ 
tions  which  are,  in  some  sense,  extremal.  Nesi  also  compares  the  finite  and  infinite  rank  schemes  and 
shows,  in  particular,  that  the  latter  is  more  powerful. 

These  results  on  attainability  strongly  suggest  that  some  bounds  on  the  effective  moduli  of 
polycrystalline  composites  remain  to  be  discovered.  Such  new  bounds  are  now  being  sought,  by  trying  to 
make  use  of  additional  null-Lagrangians  besides  those  used  in  [7], 
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(17)  Iterated  Homogenization  and  Differential  Effective  Medium  Theory 
(Avellaneda,  Milton,  Papanicolaou) 

In  1935  Bruggeman  [1]  introduced  a  class  of  dielectrics  for  which  the  effective  or  macroscopic  pro¬ 
perties  could  be  computed  exactly  from  a  relatively  simple  formula.  His  theory  came  to  be  known  as  the 
“Effective  Medium  Theory”  and  is  used  extensively  in  a  variety  of  applications.  Bruggeman  described 
the  class  of  materials  for  which  his  theory  was  valid  in  a  qualitative  way:  they  have  to  contain  inhomo¬ 
geneities  with  arbitrarily  small  scales,  and  the  addition  of  inhomogeneities  of  one  scale  can  be  dealt  with 
incrementally,  as  if  all  inhomogeneities  at  coarser  scales  were  homogenized  or  had  only  macroscopic 
effects.  This  almost  axiomatic  description  of  the  materials,  followed  by  the  formula  for  their  effective 
dielectric  constant  with  virtually  no  other  computations,  was  a  veritable  tour  de  force  by  a  brilliant  physi¬ 
cist.  But  because  the  formula  that  he  postulated  is  a  simple  one  and  because  it  can  also  be  derived  (as  was 
done  in  the  1960’s)  by  a  self-consistent  scheme,  confusion  arose  as  to  precisely  what  class  of  materials 
has  this  effective  behavior.  More  mathematically  inclined  researchers  can  go  a  step  further  and  question 
the  existence  of  real  materials  with  Bruggeman’s  effective  behavior  altogether. 

Milton  [2]  first  showed  that  Bruggeman’s  theory  (called  the  coherent  potential  approximation  in  the 
meantime)  is  realizable.  His  proof  is  complicated  and  leaves  the  reader  with  the  general  impression  that 
the  class  of  materials  for  which  the  theory  works  is  small.  Norris  [3]  gives  a  contemporary  account  of 
effective  medium  theory,  including  the  differential  effective  medi  ■  n  theory  in  which  the  volume  occu¬ 
pied  by  the  inhomogeneities  of  the  various  scales  need  not  be  all  the  material  volume.  The  conceptual 
framework  is  a  more  qualitative  version  of  Bruggeman’s  rules  or  “axioms”  but  the  issue  of  what  is 
exactly  the  underlying  material  is  not  addressed. 

Avellaneda  [4]  gives  a  very  simple  and  rather  complete  mathematical  description  of  effective 
medium  theory  starting  with  real  materials  and  using  the  theory  of  iterated  homogenization  [5]. 
Bruggeman’s  “axioms”  are  realized  nicely,  as  they  were  in  Milton’s  work  [2],  and  it  is  now  clear  that  a 
very  broad  class  of  materials,  called  hierarchical  materials,  can  be  described  by  effective  medium  theory 
and  its  gen^.'  jiizations.  Iterated  homogenization  exploits  very  efficiently  the  separation  of  scales  present 
in  the  hierarchical  materials.  The  ideas  here  are  closely  connected  with  the  renormalization  group  formal¬ 
ism  where  the  separation  of  scales  is  done  in  a  heuristic  way.  In  fact  the  renormalization  group  formalism 
applied  to  composite  dielectrics  reduces  to  iterated  homogenization  if  separation  of  scales  is  imposed. 
The  difficult,  deep  and  open  question  here  is  why  the  renormalization  group  works  even  when  there  is  no 
clear  separation  of  scales. 
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In  [6]  Avellaneda  and  Papanicolaou  analyze  a  class  of  nonlinear  composites  with  hierarchical  struc¬ 
ture  and  generalize  the  effective  medium  theory  to  fully  nonlinear  materials.  Nonlinear  iterated  homogen¬ 
ization  is  used  and  both  the  symmetric  (Bruggeman)  version  as  well  as  the  differential  effective  medium 
theory  are  studied. 


References 

[1]  D.  A.  G.  Bruggeman,  “Berechnung  verschiedener  physikalischer  Konstante  von  heterogene 
Substanzen”,  Ann.  Physik,  24,  636  (1935). 

[2]  G.  W.  Milton,  “The  coherent  potential  approximation  is  a  realizable  effective  medium 
scheme”,  Comm.  Math.  Phys.,  99, 463-500  (1985). 

[3]  A.  N.  Norris,  “A  differential  scheme  for  the  effective  moduli  of  composites”,  Mechanics  of 
Materials,  4, 1-16,(1985). 

[4]  M.  Avellaneda,  “Iterated  homogenization,  differential  effective  medium  theory  and  its  applica¬ 
tions”,  Comm.  Pure  Appl.  Math.,  40,  527-554  (1987). 

[5]  A.  Bensoussan,  J.  L.  Lions  and  G.  Papanicolaou,  Asymptotic  Analysis  of  Periodic  Structures, 
North-Holland,  1978. 

[6]  M.  Avellaneda  and  G.  Papanicolaou,  “Iterated  homogenization  and  nonlinear  effective  medium 
theories”,  submitted  to  Comm.  Pure  Appl.  Math.,  1989. 


i 


-59- 


(18)  Numerical  Approximations  of  Oscillatory  Solutions  to  Partial  Differential  Equations 
(Avellaneda,  Hou,  Papanicolaou) 

There  are  many  computational  problems  with  highly  oscillatory  solutions,  e.g.  computation  of  high 
frequency  acoustic  and  electro-magnetic  fields,  properties  of  composite  materials  and  turbulent  flow. 
Sometimes  the  original  equations  describing  this  problems  can  be  replaced  by  effective  equations  model¬ 
ing  some  average  quantity  without  the  oscillations.  Geometrical  optics,  homogenization  of  composite 
materials  and  turbulence  models  are  such  examples  [3]. 

Whenever  effective  equations  are  applicable,  they  are  also  very  useful  for  computational  purposes. 
However,  there  are  many  situations  for  which  we  do  not  have  well-posed  effective  equations  or  for  which 
the  solution  contains  different  frequencies  such  that  effective  equations  are  not  practical.  In  these  cases, 
the  question  naturally  arises  whether  the  original  formulation  can  be  used  to  obtain  useful  approximations 
without  accurately  resolving  small  scale  structures  on  the  computational  grid. 

In  this  project,  we  investigate  the  possibility  of  approximating  original  equations  directly  without 
resolving  all  scales  in  the  solutions.  The  classes  of  equations  we  consider  include  semilinear  hyperbolic 
systems,  elliptic  problems  and  incompressible  Euler  equations  with  oscillatory  solutions.  In  particular, 
we  study  convergence  properties  of  several  numerical  methods,  i.e.  finite  difference  methods,  particle 
methods  and  random  choice  methods,  when  applied  to  these  problems  with  oscillatory  solutions. 

There  are  three  main  sources  of  difficulties  for  discrete  approximations  of  highly  oscillatory  solu¬ 
tions. 

(i)  The  first  is  the  sampling  of  the  computational  mesh  points  (xj  =j  Ax ,  y'=0, 1,...).  There  is  the 
risk  of  resonance  between  the  mesh  points  and  the  oscillations.  If,  e.g.  Ax  equals  the  wave  length  of  the 
oscillation,  the  discrete  initial  data  may  only  get  values  from  the  peaks  of  the  oscillatory  solution.  We  can 
never  expect  convergence  in  that  case.  Thus  Ax  cannot  be  completely  independent  of  the  wave  length  [8- 
10], 

(ii)  Another  difficulty  comes  from  the  approximations  of  advection.  The  group  velocity  for  the  dif¬ 
ferential  equation  and  the  corresponding  discretization  are  often  very  different  This  means  that  an  oscil¬ 
latory  pulse  which  is  not  well  resolved,  is  not  transported  correctly  by  the  approximation.  Furthermore, 
dissipative  schemes  rarely  advect  oscillations  correctly.  The  oscillations  are  damped  out  very  fast  [8,1 1]. 

(iii)  Finally,  the  nonlinear  interaction  of  different  high  frequency  components  in  a  solution  must  be 
modeled  correctly.  High  frequency  interaction  may  produce  lower  frequencies  that  influence  the 
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averaged  solution  [8,13,16]. 

Bjorn  Engquist  was  the  first  to  study  this  problem.  He  showed  that  for  the  semilinear  Carleman 
equations  [5]  with  oscillatory  initial  data,  a  particle  method  converges  essentially  independent  of  the 
oscillatory  wave  length  (say  e)[8].  Convergence  essentially  independent  of  e  means  that  the  error  is  small 
for  all  samplings  of  the  grid  outside  a  set  of  arbitrarily  small  measure  [8]. 

Tom  Hou  continued  the  study  of  this  problem  in  his  Ph.D.  thesis  under  the  supervision  of  Engquist 
[14].  They  chose  the  semilinear  discrete  Boltzmann  equations  in  kinetic  theory  as  a  model  problem  [15]. 
In  particular,  they  consider  the  Broadwell  model  whose  properties  are  very  different  from  those  of  the 
Carleman  model  [4],  In  such  equations,  the  high  frequency  components  can  also  be  transformed  to  the 
lower  frequencies  through  the  nonlinear  interaction,  which  depend  sensitively  on  the  velocity  components 
[13,16].  Their  study  shows  that  particle  methods  can  provide  a  good  approximation  to  highly  oscillatory 
solutions  of  the  discrete  Boltzmann  equations.  The  crucial  step  in  their  particle  method  is  the  updating  of 
the  numerical  solutions  along  their  characteristics  without  using  local  interpolation.  This  guarantees  that 
the  advection  step  is  correctly  approximated  without  introducing  numerical  dispersion  or  dissipation. 
Moreover,  since  the  oscillations  in  the  solutions  are  propagated  along  characteristics,  the  particle  method 
approximation  captures  the  nonlinear  interactions  properly.  Thus  if  enough  ergodic  mixing  is  provided  in 
the  samplings  of  the  grid,  classical  ergodic  theory  implies  that  the  numerical  solution  will  converge  to  the 
continuous  averaged  solution  with  probability  one  [1,8,9,10].  That  is,  the  set  containing  bad  samplings 
has  arbitrarily  small  measure. 

Another  application  .o  Lhc  vertex  method  [6]  for  2-P  incompressible  Euler  equations  with  highly 
oscillatory  vorticity.  Although  many  numerical  experiments  seem  to  suggest  that  vortex  methods  could 
produce  useful  results  for  high  Reynolds  number  flow  without  resolving  the  small  scales  in  the  solution 
(see  e.g.  [17]),  there  is  still  lack  of  theoretical  justification  for  these  observations.  For  our  model  problem 
with  oscillatory  but  bounded  vorticity  field,  we  can  show  that  the  vortex  method  converges  to  the  homo¬ 
genized  solutions  of  the  Euler  equations  essentially  independent  of  the  oscillation  [7,14],  The  analysis 
here  follows  from  similar  observations  for  the  discrete  Boltzmann  equations. 

There  are  problems  with  oscillatory  solutions  for  which  sampling  of  the  grid  alone  is  not  enough  to 
guarantee  convergence  of  the  numerical  method  essentially  independent  of  the  oscillation.  For  second 
order  elliptic  equations  with  oscillatory  coefficients,  Tom  Hou  in  collaboration  with  Marco  Avellaneda 
and  George  Papanicolaou  shows  that  except  in  one  dimensional  cases,  finite  difference  methods  in  gen¬ 
eral  will  not  converge  to  the  correct  homogenized  solutions  with  probability  one  independent  of 
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oscillations  [2].  Certain  small  scale  features  in  the  solution  have  to  be  resolved  well  by  the  computational 
grid.  On  the  other  hand,  it  is  still  possible  to  obtain  useful  approximations  without  resolving  all  the  small 
scale  features  in  the  solutions.  In  the  case  of  periodic  coefficients,  we  show  that  it  is  enough  to  require 
(h/e-[h/e])  to  be  small,  where  h  is  the  grid  size,  e  is  the  oscillatory  wavelength  and  [h/c]  is  the  integer 
part  of/i/e  [2], 

In  the  near  future  we  would  like  to  extend  the  sampling  idea  to  problems  in  fluid  dynamics  where 
the  oscillations  are  produced  by  the  equations  themselves.  This  leads  to  the  well-known  subgrid  model¬ 
ing  problem  for  which  homogenization  ideas  are  only  marginally  useful  [18]  at  present. 
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(19)  Recent  Progress  in  Traveltime  Tomography  (Berryman) 

One  recent  innovation  in  traveltime  tomography  is  the  introduction  of  variational  methods  to  aid  in 
the  construction  of  the  initially  "close"  models  of  structure  required  for  the  success  of  diffraction  tomog¬ 
raphy  in  all  but  the  simplest  situations.  The  essential  idea  is  to  construct  a  variational  functional  based  on 
the  physics  of  wave  propagation  through  elastic  materials,  i.e.,  the  minimum  of  the  variational  functional 
occurs  when  the  waves  satisfy  the  correct  partial  differential  equations  [1,2].  The  tomographic  data 
obtained  from  all  source  and  receiver  pairs  then  form  a  complicated  set  of  constraints  on  the  variation 
which  are  easily  treated  on  the  same  footing  as  the  wave  equations.  The  advantage  of  this  approach  is 
that,  instead  of  solving  the  partial  differential  equations  for  the  current  (approximate)  model  exactly 
everywhere  at  each  step  of  an  iteration  process,  variational  methods  can  be  used  to  obtain  reasonable 
approximate  solutions  quickly.  For  example,  there  is  clearly  not  much  point  in  satisfying  Snell’s  law 
exactly  at  every  interface  for  a  model  structure  whose  chosen  cell  pattern  is  arbitrary  and  whose  wave 
speeds  are  initially  expected  to  be  a  poor  approximation  to  the  true  values;  in  particular,  since  the  loca¬ 
tions  of  interfaces  between  cells  are  generally  chosen  in  an  arbitrary  manner,  the  choices  made  should  not 
have  a  large  effect  on  the  final  result,  but  with  common  methods  often  do.  The  new  approach  that  has 
been  developed  by  Berryman  [3-5]  constructs  rays  that  satisfy  Fermat’s  principle  globally  and  approxi¬ 
mately  satisfy  Snell’s  law  locally  at  each  stage  of  the  calculation  and  with  a  level  of  accuracy  consistent 
with  that  of  the  model.  The  method  has  been  implemented  in  two-dimensions  by  Berryman  and  Yorkey 
[6]  using  a  new,  simple,  and  fast  simplex  algorithm  [7]  for  finding  optimum  coefficients  in  the  Fourier 
sine  series  expansion  of  the  graph  of  the  ray.  This  ray  tracing  method  is  currently  in  use  in  our  traveltime 
reconstructions,  and  provides  an  efficient  means  of  generating  the  "close"  model  of  structure  required  as 
input  to  the  sophisticated  methods  of  Devaney,  which  are  expected  to  converge  more  rapidly  once  they 
can  be  applied.  In  some  applications,  the  model  obtained  from  traveltime  tomography  may  provide  a 
satisfactory  estimate  of  the  wave  speed  structure  by  itself,  without  needing  the  additional  data  and  effort 
required  for  diffraction  tomography. 

New  algorithms  have  recently  been  developed  for  inverting  traveltime  data  in  media  with  high 
wave  speed  contrasts  (and  therefore  with  strong  ray  bending)  [3,4],  One  of  the  most  significant  new  con¬ 
tributions  of  this  work  has  been  the  realization  that  Fermat’s  principle  (i.e.,  the  true  ray  path  for  the  first 
arrival  is  the  one  of  least  traveltime)  makes  a  very  powerful  statement  about  the  feasibility  or  lack  of 
feasibility  of  trial  rays  through  a  given  wave  speed  model  for  the  medium  of  interest  [4-6].  Nonfeasible 
rays  are  trial  rays  with  traveltimes  shorter  than  that  measured  —  such  rays  are  physically  impossible. 
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Feasible  rays  are  all  other  rays  with  traveltimes  as  long  as  or  longer  than  that  measured.  Fermat’s  princi¬ 
ple  may  therefore  be  used  to  "convexify"  the  fully  nonlinear  inversion  problem.  Since  there  exists  a 
definite  convex  region  of  the  slowness  model  space  —  depending  only  on  the  measured  traveltimes  -  that 
contains  all  feasible  models,  the  "solution"  (if  one  exists)  to  the  traveltime  inversion  problem  lies  on  the 
boundary  of  this  convex  set.  Thus,  traveltime  tomography  can  be  reduced  to  a  problem  in  nonlinear  con¬ 
strained  optimization.  Obtaining  an  exact  solution  of  this  optimization  problem  would  be  very  hard 
indeed  because  the  pertinent  vector  space  is  multidimensional  and,  although  every  point  on  the  feasibility 
surface  can  be  found  in  principle,  the  computations  involved  in  mapping  such  a  multidimensional  surface 
are  prohibitive.  However,  our  knowledge  of  the  existence  of  this  feasibility  surface  is  sufficient  to  pro¬ 
vide  a  simple,  approximate  figure  of  merit  that  allows  us  to  stabilize  the  reconstruction  algorithm  with  a 
negligible  increase  in  the  computation  time. 


References 

[1]  J.R. Willis,  “Variational  principles  for  dynamic  problems  for  inhomogeneous  elastic  media,” 
Wave  Motion  3, 1-11  (1981). 

[2]  M.E.Gurtin,  “Variational  principles  for  linear  elastodynamics,”  Arch.  Rational  Mech.  Anal.  16, 
34-50  (1964). 

[3]  J.  G.  Berryman,  “Weighted  least-squares  criteria  for  seismic  traveltime  tomography,”  IEEE 
Trans.Geosci  .Remote  Sensing,  to  appear  May,  1989. 

[4]  J.  G.  Berryman,  “Traveltime  tomography  and  nonlinear  constrained  optimization,”  to  appear  in 
Proceedings  of  the  Twenty-Second  Annual  Asilomar  Conference  on  Signals,  Systems,  and  Com¬ 
puters,  Oct.31-Nov.2, 1988,  Pacific  Grove,  CA 

[5]  J.  G.  Berryman,  “New  applications  of  Fermat’s  principle  to  seismic  traveltime  tomography,” 
invited  talk  at  the  17th  International  Conference  on  Mathematical  Geophysics  in  Blanes,  Spain, 
June  20-25, 1988. 

[6]  J.  G.  Berryman  and  T.  J.  Yorkey,  “New  applications  of  Fermat’s  principle  to  traveltime  tomogra¬ 
phy,”  Phys.  Rev.  Lett.,  submitted.. 

[7]  J.  A.  Nelder  and  R.  Mead,  “A  simplex  method  for  function  minimization,”  Computer  J.  7, 
308-313(1965). 


-64- 


Vn.4  Nonlinear  Optics 

(20)  The  Focusing  Singularity  of  the  Nonlinear  Schroedinger  Equation 
(Landman,  LeMesurier,  Papanicolaou,  Sulem,  Sulem) 

The  nonlinear  Scroedinger  equation  (NLS) 

+  Avy  +  |v|/|2y  =  0 
\|/(0,x)  =  vjt0(x) ,  xe.RN 

arises  frequently  in  the  description  of  wave  phenomena  as  envelope  equation  for  field  amplitudes.  In  non¬ 
linear  optics  it  describes  approximately  the  propagation  of  narrow  beams  in  a  material  whose  index  of 
refraction  increases  or  decreases  with  field  intensity  leading  to  focusing  or  defocusing  of  the  beam.  Here 
the  variable  t  represents  distance  along  the  direction  of  the  beam  axis  and  the  spatial  variables 
x  =  (*!  ,x2)  e  R2  are  coordinates  in  the  cross  sectional  plane  of  the  beam.  In  plasma  physics  it  arises  with 
N  =  3  and  it  describes  the  evolution  of  Langmuir  waves  in  plasmas,  in  the  subsonic  limit.  Recent  surveys 
of  applications  and  basic  properties  of  the  nonlinear  Schroedinger  equation  are  [1,2]. 

The  main  mathematical  problem  associated  with  NLS  is  the  existence  of  singular  solutions  in  two 
or  more  dimensions  (N£ 2),  their  local  form  and  the  nature  of  the  solutions  (if  they  exist)  past  the  singular¬ 
ity.  The  singularity  occurs  when  the  nonlinear  focusing  overcomes  wave  dispersion,  the  beam  amplitude 
becomes  infinite  on  its  axis  and  its  cross  section  shrinks  to  zero.  These  basic  facts  have  been  known  for 
over  twenty  years  [1,2]  but  the  local  form  of  the  singularity  is  difficult  to  determine,  analytically  or 
numerically,  and  several  attempts  in  the  literature  proved  later  to  be  contradictory.  Understanding  of  the 
local  form  of  the  singular  solutions  is  important  for  many  reasons.  Mathematically  it  is  a  problem  com¬ 
parable  to  understanding  the  structure  of  discontinuous  solutions  of  model  conservation  laws  such  as 
Burger’s  equation  or  the  Korteweg-De  Vries  equation  in  the  zero  dispersion  limit:  it  is  a  canonical  situa¬ 
tion  in  which  competing  effects  (focusing  nonlinearity  and  wave  dispersion)  produce  phenomena  that 
appear  quite  generally  in  much  more  complex  systems.  Physically  and  numerically  it  is  an  equation  on 
which,  because  of  its  relative  simplicity,  one  can  develop  and  study  expansions,  computational  schemes, 
etc.  that  can  be  used  in  more  elaborate  settings:  saturated  nonlinearities  (where  jq/j2  is  replaced  by 
|\|/|2(1  +  |v|2)-1)  and  effects  of  backscattering  in  beams,  the  Zaharov  equations  in  plasmas  and  so  on. 
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George  Papanicolaou  in  collaboration  with  D.  McLaughlin  (University  of  Arizona)  begun  a  study 
of  the  local  form  of  the  focusing  singularity  in  1980.  Little  progress  was  made  until  1984  when  in  colla¬ 
boration  with  B.  LeMesurier  (then  a  graduate  student,  now  at  the  University  of  Arizona),  C.  Sulem  and 
P.L.  Sulem  (from  France  and  Israel)  extensive  numerical  computations  were  conducted.  Many  numerical 
computations  of  singular  solutions  of  NLS  had  of  course  been  done  by  other  investigators  previously,  but 
they  did  not  have  the  resolution  near  the  singularity  to  distinguish  between  competing  theoretical  con¬ 
structions  based  on  asymptotic  analysis.  The  numerical  computation  of  singular  solutions  of  an  equation 
in  which  there  are  no  dissipative  effects  whatsoever,  such  as  NLS  or  a  wave  equation,  is  exceedingly 
difficult.  In  [3]  we  introduced  a  method,  we  call  it  dynamic  rescaling,  which  allowed  us  to  compute 
numerically  singular  solutions  up  to  amplification  of  order  1012.  The  idea  is  to  let  the  size  of  the  numeri¬ 
cal  mesh  depend  on  time  and  to  derive  an  equation  for  it  that  has  to  be  solved  along  with  NLS  (suitably 
rescaled).  An  equation  for  the  mesh-scale  factor  is  obtained  by  demanding  that  a  suitable  norm  of  the 
wave  function  stays  bounded.  The  new  coupled  system  (NLS  plus  the  mesh-scale  factor  equation)  has  no 
singularities  and  therefore  can  be  solved  by  standard  numerical  methods.  We  used  two  different  numeri¬ 
cal  methods  to  check  our  results.  The  nature  of  the  focusing  singularity  can  be  read  off  from  the  time  evo¬ 
lution  of  the  mesh-scale  factor  which  is  computed  very  accurately. 

In  [4,5]  we  presented  more  extensive  calculations.  We  concluded  that  in  dimensions  N> 2,  the 
supercritical  case,  solutions  blow  up  like  (r*  -  where  t*  is  the  time  where  the  solution  becomes 
singular,  and  we  characterized  the  spatial  profile  of  singular  solutions  through  a  nonlinear  eigenvalue 
problem.  We  have  not  succeeded  up  to  now  in  providing  a  full  theory  for  this  difficult  mathematical  prob¬ 
lem,  but  we  have  obtained  a  lot  of  information  via  formal  asymptotics  and  numerical  computations.  In 
two  dimensions  N  =  2,  the  physically  interesting  case  in  nonlinear  optics,  our  computations  showed  that 
the  blow-up  rate  (r*-r)~2'3,  predicted  by  Zaharov  and  Synach  [8],  or  [(r*  -r)/log(r*  -f)-1r1/2. 
predicted  by  Talanov  [9]  and  Wood  [10],  are  not  correct.  We  were  not  able  however,  even  with  our  com¬ 
putations,  to  guess  the  true  form  of  the  blow-up  rate.  In  [4,5]  we  included  many  rather  extensive  calcula¬ 
tions  of  solutions  of  NLS  with  saturated  nonlinearities.  They  have  many  intriguing  features,  such  as 
dispersive  damping  after  repeated  focusings,  that  we  do  not  understand  analytically.  These  calculations 
were  done  on  a  Cray  XM-P.  The  calculation  of  the  singular  solutions  without  saturation  was  done  on  our 
own  mini  computers  (mostly  the  Celerity  1200). 


By  the  summer  of  1987  we  were  convinced  that  no  amount  of  numerical  computation,  no  matter 
how  elaborate  or  accurate,  could  help  us  determine  the  form  of  the  focusing  singularity  in  two  dimensions 


N  =  2,  the  critical  case.  Using  a  matched  asymptotic  expansion  we  derived  in  [51  a  solution  that  was 
singular  in  the  same  way  Zaharov  and  Synach  [8]  predicted,  by  very  different  methods.  In  reworking  our 
analysis  we  realized  that  our  asymptotic  solutions  (and  those  of  [8])  had  a  very  different  profile  than  the 
one  that  the  numerical  computations  produced.  There  was  no  way  to  fix  this  except  by  perturbing  the 
dimension  to  JV>2  and  then  using  the  profiles  from  the  nonlinear  eigenvalue  problem  of  the  supercritical 
case.  Those  profiles  were  much  closer  to  the  numerically  observed  profiles  when  N> 2  was  near  2.  With 
this  insight,  and  in  collaboration  with  M.  Landman,  we  were  able  to  determine  completely  [6,7]  the 
singular  blow-up  in  the  critical  case  which  is  [(/*  -  r)/log  log(f*  -  r)-1  ]-5/2.  This  also  explains  why  it  is 
difficuh  to  get  such  information  from  numerical  calculations. 

We  are  currently  working  on  an  extension  of  the  dynamic  rescaling  ideas  to  more  general  situations 
where  there  are  multiple  peaks  in  the  solution,  the  peaks  move  and  they  need  not  have  (locally)  isotropic 
profiles.  Our  aim  is  to  write  a  general  code  in  three  space  dimensions  in  which  we  implement  these  ideas. 
We  are  also  studying  intensively,  with  X.  Wang  (a  graduate  student),  the  mathematical  questions  that 
arise  in  connection  with  the  nonlinear  eigenvalue  probler  for  the  singular  profiles. 
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(21)  Nonlinearity  and  Localization  in  One  Dimensional  Random  Media 
(Knapp,  Papanicolaou,  White) 

R.  Knapp,  G.  Papanicolaou  and  B.  White  [1]  have  analyzed  the  combined  effects  of  nonlinearity 
and  randomness  in  the  propagation  of  time  harmonic  waves  in  one  dimensional  media.  In  nonlinear  opti¬ 
cal  media  the  index  of  refraction  depends  on  the  intensity  of  light,  giving  rise  to  many  interesting  effects, 
including  optical  bistability.  When  an  optical  device  can  have  two  different  output  states  for  a  given 
input  intensity  (depending  on  hysteresis)  it  is  said  to  be  bistable.  On  the  other  hand  in  linear  media  with 
random  inhomogeneities,  different  phenomena  arise,  such  as  localization.  When  localization  occurs  the 
transmitted  intensity  decays  exponentially  as  a  function  of  the  length  of  material  the  light  passes  through. 

Optical  bistabilty  was  first  observed  experimentally  in  a  Fabry-Perot  etalon.  The  etalon  is  con¬ 
structed  from  a  slab  of  nonlinear  material  of  thickness  L  sandwiched  between  two  partially  reflecting 
dielectric  mirrors.  Randomness  can  be  introduced  either  from  impurities  or  by  constructing  the  medium 
with  alternating  films  of  nonlinear  and  linear  material  of  random  thickness.  We  consider  a  time  harmonic 
model  for  the  etalon  with  a  wave  incident  from  the  right 

«xr+*2(l  +eit(x)  +  w  !u  l2)u=0  0 <x<L 
u  (x)=€_lfa  +  Re**1  x  >L 
u(x)=Te~,kx  x  <0 

The  field  amplitude  has  been  scaled  by  the  amplitude  of  the  incident  wave  so  that  the  nonlinear 
coefficient  w  is  proportional  to  the  input  intensity.  The  transmission  and  reflection  coefficients  are 
respectively  T  and  R.  The  randomness  is  modeled  by  p(z)  which  is  a  bounded,  stationary,  zero  mean  sto¬ 
chastic  process.  The  parameter  e  allows  us  to  control  the  size  of  the  noise.  Requiring  continuity  of  u  and 
uz  at  x=0  and  x=L  gives  us  a  boundary  value  problem.  We  are  considering  the  case  where  the  zero- 
intensity  homogeneous  medium  is  matched  to  free  space  (index  of  refraction  equal  to  one)  to  simplify  the 
analysis. 

A  significant  difficulty  in  the  analysis  is  the  possibility  of  bistability  or  lack  of  uniqueness.  For  a 
given  value  of  w  there  may  be  several  solutions  to  the  time  harmonic  problem.  Because  of  the  non¬ 
uniqueness  it  is  desirable  to  formulate  the  problem  in  a  slightly  different  way.  We  convert  it  into  a  family 
of  initial  value  problems  with  a  parameter  proportional  to  the  transmitted  intensity  This  allows  us  to  use 
methods  from  the  asymptotic  theory  of  stochastic  equations  when  suitable  scalings  t r;  introduced.  When 
there  is  no  randomness  we  find  that  bistability  occurs  when  h’=L-i/3  for  L  — >o°.  When  random 
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inhomogeneities  are  present  the  analysis  of  the  transmission  problem  is  much  more  difficult  and  we  first 
study  the  growth  of  the  solutions  of  the  family  of  initial  value  problems  when  the  noise  intensity  e  is 
small  and  L  is  large.  In  the  linear  case  it  is  exponential  but  in  the  nonlinear  case  it  is  only  algebraic.  This 
has  important  implications  for  the  onset  of  bistability  when  random  inhomogeneities  are  present.  It  occurs 
for  much  smaller  L,  compared  to  he  nonrandom  case.  We  have  made  a  full  study  of  the  statistical  proper¬ 
ties  of  the  onset  of  bistability  both  analytically  with  asymptotic  methods  and  numerically. 
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Dec.  1-5,  1986,  IMA  Univ.  of  Minn.,  workshop  on  math,  methods  in  exploration  geophysics,  invited 
speaker. 

Jan.  21-23,  1987,  AMS  annual  meeting,  special  session  on  homogenization  and  oscillations,  San  Antonio, 
Texas,  invited  speaker. 

Feb.  18-21, 1987,  University  of  Arizona  conference  on  random  media,  invited  speaker. 

March  4-5,  1987,  two  invited  lectures  at  Purdue  University. 

March  16-20, 1987,  RPI  conference  on  multiphase  flow,  invited  participant. 

April  1-4, 1987,  U.  of  Cambridge,  England,  three  lectures  at  the  annual  meeting  of  British  Probabilists. 
May  4,  1987,  invited  lecture  at  the  B.  Levich  Memorial  Conference  at  the  City  University  of  New  York. 
Aug.  24-28,  1987,  Equadiff,  Xanthi,  Greece,  co-organizer  and  speaker. 

Aug  24-28,  1987,  special  invited  speaker  at  the  European  Statisticians  Conference,  Thessaloniki,  Greece. 

Sept.  14-16, 1987,  Duke  University  and  University  of  North  Carolina  at  Chapel  Hill,  visit  and  lecture. 

Nov.  14-15,  1987,  invited  speaker  at  the  annual  meeting  of  the  Soc.  for  Natural  Philosophy  at  Yale  U. 

March  14-18, 1988,  Tuscon  Arizona,  meeting  on  singularities  in  PDE 

April  15,16, 1988,  Brown  University,  meeting  in  honor  of  Wendell  Reming 

May  4-6,  1988,  Los  Alamos,  Meeting  on  Nonlinearity  and  Localization 

May  18-26,  1988,  L’Aquilla  Italy,  Meeting  on  Interacting  particle  systems 

June  6-10,  1988,  Paris  France,  Meeting  in  honor  of  J.  L.  Lions 

July  18-29,  1988,  Silivri  Turkey,  Workshop  on  stochastic  analysis. 

Oct  22-25,  1988,  Cornell  University,  Workshop  on  multiphase  flow. 

Dec  6-9, 1988,  Workshop  on  composite  and  random  media,  Leesburg  Va. 

Dec  14-16, 1988,  Conference  on  computational  methods  in  wave  propagation,  Nice,  France. 


Home  Address: 


110  Bleecker  Street.  #14-CS 
New  York,  NY  10012 
(212)  533-0390  (home) 
(212)  998-3139  (office) 


Education: 

Chiao  Tung  University,  BSc,  1946 
MIT,  MSc  1948 
Harvard  University,  MSc,  1949 
New  York  University,  DSc,  1951 

Employment: 

1951- 52,  Research  Associate,  New  York  University 

1952- 55,  Special  Design  Engineer,  Foster  Wheeler  Corp. 
1955-57,  Research  Associate 

1957- 58,  Research  Assistant  Professor 

1958- 60,  Research  Associate  Professor 

1960-64,  Research  Professor  Polytechnic  Institute  of  Brooklyn 
1964-68,  Professor  of  Aeronautics  and  Astronautics 
1968-Present,  New  York  University 
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Memberships: 

Sigma  Xi 

Sigma  Gamma  Tau 

Society  for  Industrial  and  Applied  Mathematics 
American  Physical  Society 
American  Institute  of  Aeronautics  and  Astronautics 
New  York  Academy  of  Sciences 


Recent  Invited  Talks  of  Lu  Ting 

1.  Minisymposium/28  "  Analytical  and  Numerical  Solutions  of  a  Vortical  Flow",  Chairperson  and 
Organizer  L.  Ting,  in  the  First  Intern.  Conf.  on  Industrial  and  Appl.  Math.,  Paris,  July  1, 1987. 

2.  L.  Ting,  "Appropriate  Boundary  Conditions  for  Problems  in  Wave  Propagation  and  in  Fluid 
Dynamics"  invited  seminar.  May  15,  1987  ,Air  Force  Wright  Aeronautical  Laboratory,  Dayton,  Ohio. 

3.  L.  Ting,  "  Dynamics  of  Three  Vortices,  Revisited  "  Applied  Mathematics  Seminar,  Stanford 
University,  Dec.  13,  1987. 

4.  L.  Ting,  10  Lectures  series  "Motion  and  Interaction  of  Vortices  with  Diffusive  Core  Structures" 
in  commemoration  of  75th  year  of  the  Aerodynamisches  Institut,  RWTH  Aachen,  April  25-29, 1988. 

5.  M.  J.  Miksis  and  L.  Ting  "  Effects  of  Bubbly  Layers  on  Wave  Propagation",  A  Workshop  on 
Theoretical  Aspects  of  Multiphase  Flow  Phenomena,  Cornell  University  Mathematical  Sciences  Institute, 
Oct.  22-25,  1988. 

6.  L.  Ting  and  M.  J.  Miksis,  "  On  Vortical  Flow  and  Sound  Generation",  EUROMECH  Collo¬ 
quium  247,  Vortex  Interaction  Noise,  Goettingen,  June  28  -  July  1, 1989. 

7.  L.  Y.  Cheng  and  L.  Ting, "  Void  Distribution  in  a  Shallow  Bubbling  Pool",  26th  National  Heat 
Transfer  Conference,  Philadelphia,  August  6-9,  1989. 
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IX.  Publications  of  the  Center 


James  G.  Berryman 

[1]  James  G.  Berryman,  Graeme  W.  Milton  -  "Microgeometry  of  Random  Composites  and  Porous 
Media"  J.  Phys.  D:  Appl.  Phys.  21  (1988)  87-94. 

[2]  James  G.  Berryman,  Lawrence  Livermore  National  Laboratory  -  "Weighted  Least-Squares  Criteria 
for  Seismic  Traveltime  Tomography"  -  preprint  Geophys.  J.  Roy.  Astron.  Soc. 


Russel  E.  Caflisch 

[1]  Lectures  on  the  Mathematical  Theory  of  Multi-Phase  Flow  (with  J.  Rubinstein).  Cl  MS  Lecture 
Note  Series  (1986). 

[2]  "Equilibrium  for  Radiation  in  a  Homogeneous  Plasma"  (with  C.D.  Levermore),  Phys.  Fluids,  29 
(1986)  748-752. 

[3]  "Mathematical  Method  for  the  Equations  of  Kinetics  of  Gases,"  (with  C.V.M.  Van  der  Mee)  Tran¬ 
sport  Th.  Stat.  Physics  15  (1986)  353-367. 

[4]  "Stability  of  Shock  Wave  Profiles  for  a  Model  Boltzmann,"  Equation  (with  T.P.  Liu)  15th  RGD 
Symposium,  Grado,  (1986)  3-9. 

[5]  "Statistical  Properties  of  Sedimenting  Suspensions  and  the  Point  Particle  Approximation,"  (with 
J.H.C.  Luke)  in  Advances  in  Multiphase  Flow  and  Related  Problems,  ed.  Papanicolaou,  SIAM 
(1986)  23-34. 

[6]  "The  Milne  and  Kramers  Problems  for  the  Boltzmann  Equation  of  a  Hard  Sphere  Gas"  (with  C. 
Bardos,  B.  Nicolaenko),  Comm.  Pure  Appl.  Math.,  39  (1986)  323-352. 

[7]  "Long  Time  Existence  for  Slightly  Perturbed  Vortex  Sheets"  (with  O.  Orellana),  Comm.  Pure  Appl. 
Math.,  39  (1986)  807-838. 

[8]  "Asymptotic  Expansions  of  Solutions  for  the  Boltzmann  Equation,"  J.  Transport  Th.  Slat.  Phys.,  16 
(1987)  701-725. 

[9]  "Statistical  Modeling  of  a  Dilute  Fluid-Particle  Suspension,"  in  Hydrodynamic  Behavior  and 
Interacting  Particle  Systems,  IMA  vol.  9,  ed.  Papanicolaou,  Springer  (1987)  1-12. 
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[10]  "Asymptotics  of  the  Boltzmann  Equation  and  Fluid  Dynamics,"  in  Kinetic  Theory  and  Gas  Dynam¬ 
ics  CISM  No.  239  (1988)  ed.  C.  Cercignani,  pp.  95-133. 

[11]  "Long  Time  Existence  and  Singularity  Formation  for  Vortex  Sheets,"  Workshop  on  Vortex 
Methods,  ed.  Anderson  and  Greengard. 

[12]  "Nonlinear  Analysis  for  the  Evolution  of  Vortex  Sheets,"  in  Vortex  Motion,  Proc.  IUTAM  Sympo¬ 
sium  on  Fundamental  Aspects  of  Vortex  Motion,  ed.  H.  Hasimoto  and  T.  Kambe,  North-Holland 
(1988)  75-77  and  in  Fluid  Dyn.  Research  3  (1988)  75-77. 

[13]  "Well-Posedness  of  Initial  and  Boundary  Value  Problems  for  the  Boltzmann  Equation,”  Proc.  16th 
Intern.  Symp.  Rarefied  Gas  Dyn.  (1988) 

[14]  "Singularity  Formation  for  Vortex  Sheets  and  Hyperbolic  Equations,"  3rd  Intern.  Wikshp.  on  Math. 
Aspects  Fluids  and  Plasmas  (1988) 

[15]  "Mathematical  Analysis  for  Vortex  Dynamics,"  Proc.  Workshop  on  Math.  Aspects  of  Vortex  Dyn. 
(1988)  SIAM. 

[16]  "Stability  of  Shock  Waves  for  the  Broadwell  Equations"  (with  T.-P.  Liu)  Comm.  Math.  Phys.  1J4 
(1988)  103-130. 

[17]  "Periodic  Solutions  for  Three  Sedimenting  Spheres"  (with  J.  Luke,  C.  Lim  and  A.  Sangani)  Phys. 
FI.  31  (1988)  3175-3179. 

[18]  "Singularity  Formation  and  Hl-Posedness  for  Vortex  Sheets"  (with  O.  Orellana)  SIAM  J.  Math. 
Anal,  to  appear 

[19]  "Convergence  of  the  Vortex  Method  for  Vortex  Sheets"  (with  J.  Lowengrub)  SIAM  J.  Num.  Anal, 
to  appear. 

[20]  "A  Localized  Approximation  for  Vortical  Flows"  (with  O.  Orellana  and  M.  Siegel)  SIAM  J. 
Applied  Math.,  submitted. 

[21]  "Numerical  Study  of  a  Viscous  Conservation  Law  of  Mixed  Type"  (with  M.  Affouf)  SIAM  J. 
Applied  Math.,  submitted. 

[22]  Proceedings  of  Workshop  on  Mathematical  Aspects  of  Vortex  Dynamics,  (editor)  SIAM  (1988). 
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T.  Y.  Hou 

[1]  T.  Y.  Hou,  "Convergence  of  Particle  Methods  for  Euler  and  discrete  Boltzmann  Equations  with 
Oscillatory  Solutions,"  Ph.  D.  Thesis,  UCLA,  1987. 

[2]  T.  Y.  Hou,  "  Homogenization  for  Semilinear  Hyperbolic  Systems  with  Oscillatory  Data,"  Comm. 
Pure  and  Applied  Math.,  Vol.  XLI 471-495  (1988). 

[3]  B.  Engquist  and  T.  Y.  Hou,  "  Particle  Method  Approximation  of  Oscillatory  Solutions  to  Hyper¬ 
bolic  Differential  Equations,"  accepted  by  SIAM  J.  Numerical  Analysis. 

[4]  B.  Engquist  and  T.  Y.  Hou,  "  Computation  of  Oscillatory  Solutions  to  Hyperbolic  Equations  Using 
Particle  Methods,"  Lecture  Notes  in  Mathematics  #1360,  Eds.  C.  Anderson  and  C.  Greengard, 
Springer-Verlag,  1987,  pp.  68-82. 

[5]  J.  Goodman,  T.  Y.  Hou  and  J.  Lowengrub,  "  The  Convergence  of  the  Point  Vortex  Method  without 
Smoothing,"  preprint. 

[6]  T.  Y.  Hou  and  J.  Lowengrub,  "  The  Convergence  of  the  Point  Vortex  Method  for  3-D  Euler’s  Equa¬ 
tions,"  preprint 

[7]  W.  E  and  T.  Y.  Hou,  "  Homogenization  and  Convergence  of  the  Vortex  Method  for  2-D  Euler 
Equations  with  Oscillatory  Vorticity  Fields,"  preprint,  submitted  to  Math.  Comp.. 

[8]  T.  Chacon-Rebello  and  T.  Y.  Hou,  "  A  Lagrangian  Finite  Element  Method  for  2-D  Euler’s  Equa¬ 
tions,"  preprint,  submitted  to  Comm.  Pi’ re  and  Applied  Mathematics. 

[9]  T.  Y.  Hou,  "  Convergence  of  Variable  Blob  Vortex  Methods  for  Euler  and  Navier-Stokes  Equa¬ 
tions,"  preprint 

[10]  M.  Avellaneda,  T.  Y.  Hou  and  G.  Papanicolaou,  "  Numerical  Homogenization  for  Second  Order 
Elliptic  Equations  with  Highly  Oscillating  Coefficients,"  preprint. 

[11]  T.  F.  Chan,  T.  Y.  Hou  and  P.  L.  Lions,  "  Geometry  Independent  Convergence  Results  for  Domain 
Decomposition  Algorithms,"  preprint. 

[12]  T.  F.  Chan  and  T.  Y.  Hou,  "  Domain  Decomposition  Preconditioners  for  General  Second  Order 
Elliptic  Problems,"  CAM  Report  88-16,  submitted  to  SIAM  J.  Sci.  and  StaL  Comp.. 

[13]  J.  Goodman,  T.  Y.  Hou  and  J.  Lowengrub,  "  The  Convergence  of  Vortex  Blob  Methods  with  Non¬ 
overlapping  Blobs  for  the  Euler  Equations,"  in  preparation. 
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[14]  T.  Y.  Hou  and  P.  Lax,  "  Dispersive  Approximations  in  Fluid  Dynamics,"  in  preparation. 

[15]  R.  Caflisch,  N.  Ercolani  and  T.  Y.  Hou,  "  Singularity  Formation  in  Vortex  Sheets,"  in  preparation. 

[16]  T.  Hou,  D.  Papageorgiou,  G.  Papanicolaou  and  J.  Zhu,  "  Numerical  Solution  of  2-D  Boundary 
Layer  Problems  Using  Variable  Blob  Vortex  Methods,"  in  preparation. 

[17]  T.  Chacon-Rebello  and  T.  Y.  Hou,  in  preparation. 

[18]  J.  Goodman  and  T.  Y.  Hou,  in  preparation. 

Robert  K  rasny 

[1]  "Numerical  Simulation  of  Vortex  Sheet  Evolution,"  Fluid  Dynamics  Research  3  (1988)  93-97. 

Michael  J.  Landman 

[1]  S.  Childress,  M.J.  Landman  and  H.  Strauss  (1988),  "Analysis  of  Euler  and  Navicr-Stokes  Rows 
with  Helical  Symmetry”  (pre-print). 

[2]  M.J.  Landman  (1989),  "Time  Dependent  Helical  Waves  in  Rotating  Pipe  Flow"  (to  be  submitted  to 
the  Journal  of  Fluid  Mechanics). 

[3]  M.J.  Landman  (1989),  "On  the  Search  for  Helical  Waves  in  Pipe  Poiseuille  Flow"  (in  preparation). 


Graeme  W.  Milton 

[1]  "Multicomponent  composites,  impedance  networks,  and  new  types  of  continued  fraction  1",  Comm. 
Math.  Phys.  Ill,  281-327  (198'7). 

[2]  "Multicomponent  composites,  impedance  networks,  and  new  types  of  continued  fraction  II", 
Comm.  Math.  Phys.  Ill,  329-372  (1987). 

[3]  Optimal  bounds  for  conduction  in  two-dimensional,  multiphase,  polycrystalline  media",  (with  G. 
Francfort),  J.  Stat.  Phys.  46,  161-177  (1987). 

[4]  "Transport  properties  of  touching  cylinder  pairs  and  of  the  square  array  of  touching  cylinders", 
(with  R.C.  McPhedran)  Proc.  Roy.  Soc.  Lond.,  A  411,  313-326  (1987). 

[5]  "Microgeometry  of  random  composites  and  porous  media",  (with  J.G.  Berryman),  J.  Phys.  D:  Appl. 
Phys.  21,87-94  (1988). 
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[6]  "Asymptotic  studies  of  closely  spaced,  highly  conducting  cylinders",  (with  R.C.  McPhedran  and  L. 
Poladian)  Proc.  Roy.  Soc.  Lond.  A  415,  185-196  (1988). 

[7 ]  "On  the  effective  conductivity  of  polycrystals  and  a  three-dimensional  phase-interchange  inequal¬ 
ity”  (with  M.  Avellaneda,  A.V.  Cherkaev,  and  K.A.  Lurie)  J.  Appl.  Phys.  63, 4989-5003  (1988). 

[8]  "Variational  bounds  on  the  effective  moduli  of  anisotropic  composites"  (with  R.V.  Kohn)  J.  Mech. 
Phys.  Solids,  to  appear. 

[9]  "Optimal  bounds  on  the  effective  bulk  modulus  of  polycrystals  (with  M.  Avellaneda)  SIAM  J. 
Appl.  Math,  to  appear. 

[10]  "The  classical  Hall  effect  in  two-dimensional  composites:  a  characterization  of  the  set  of  realizable 
effective  conductivity  tensors"  Phys.  Rev.  B.  to  appear. 

[11]  "On  the  conductivity  of  polycrystals  and  a  phase-interchange  inequality"  (with  M.  Avellaneda,  K. 
Lurie  and  A.  Cherkaev),  to  appear  in  Proc.  of  ETOPIM  2,  Paris,  Physica  A,  February  1989. 

[12]  "Bounds  on  the  effective  elasticity  tensor  of  composites  based  on  two-point  correlations"  (with  M. 
Avellaneda)  submitted  to  ASME  Conference  on  Composite  Materials,  Houston,  January  1989. 

[13]  "The  theory  of  field  mechanics  and  its  application  to  characterizing  the  set  of  possible  effective  ten¬ 
sors  of  composites",  submitted  to  Comm.  Pure  Appl.  Math. 


George  C.  Papanicolaou 

[1]  "Statistics  for  pulse  reflection  from  a  randomly  layered  medium"  (with  R.  Burridge  and  B. 
White)  SIAM  J.  Appl.  Math.  47  (1987)  pp.  146-168.  (Short  version  in  Lecture  Notes  in  Math. 
1186,  L.  Arnold  editor,  Springer,  1986  pp.  265-270.) 

[2]  "Asymptotic  analysis  of  the  Lyapounov  exponent  and  rotation  number  of  the  random  oscillator 
and  applications”  (with  Arnold  and  Wihstutz)  SIAM  J.  Appl.  Math.  46  (1986)  pp.  427-450 

[3]  "Focusing  singularity  of  the  nonlinear  Schroedinger  ecuation"  (with  D.  McLaughlin,  P.-L.  Sulem 
and  C.  Sulem)  Phys.  Rev.  A  34  Number  2  (1986)  pp.  1200-1 2 10. 

[4]  "Self  diffusion  of  interacting  Brownian  particles"  (with  M.  Z.  Guo)  in  "Probabilistic  Methods  in 
Mathematical  Physics"  edited  by  K.  Ito  and  N.  Ikcda,  Kinokuniya  ,  Tokyo  1987,  113-151. 
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[5]  "The  spectnim  of  backscatter  for  pulse  reflection  from  a  randomly  stratified  medium"  (with  R. 
Burridge  and  B.  White)  Springer  Lecture  Notes  IMA  (Univ.  Minnesota)  vol.  7  (1987)  edited  by 
G.  Papanicolaou,  pp.  29^14. 

[6]  "One  dimensional  wave  propagation  in  a  highly  discontinuous  medium"  (with  R.  Burridge  and  B. 
White)  Wave  Motion,  10  (1988)  pp.  19-44. 

[7]  "Multiple  scattering  noise  in  one  dimension:  universality  through  localization  length  scaling" 
(with  P.  Sheng,  B.  White  and  Z.-Q.  Zhang),  Phys.  Rev.  Letters  57  (1986)  pp.  1000-1003. 

[8]  "Minimum  wave  localization  length  in  a  one  dimensional  random  medium"  (with  P.  Sheng,  B. 
White  and  Z.-Q.  Zhang)  Phys.  Rev.  B  34  #7  (1986)  pp.  4757-4761 

[9]  "The  focusing  singularity  of  the  nonlinear  Schroedinger  equation"  (with  B.LeMesurier, 
C.Sulem  and  P.-L.  Sulem),  in  Directions  in  Partial  Differential  Equations,  Academic  Press,  1987, 
159-201.  M.  G.  Crandall,  P.  H.  Rabinowitz  and  R.  E.  L.  Turner  editors. 

[10]  "The  point  interaction  approximation  for  diffusion  in  regions  with  many  small  holes"  (with 
R.Figari  and  J.Rubinstein)  Lecture  Notes  in  Biomathematics,  volume  70,  Springer,  1987,  202-209. 
M.  Kimura,  G.  Kallianpur  and  T.  Hida  editors. 

[11]  "Direct  and  inverse  problems  for  pulse  reflection  from  a  stratified  random  halfspace  "  (with  R. 
Burridge,  P.  Sheng  and  B.  White)  in  Advances  in  multiphase  flow  and  related  problems,  G.  Papan¬ 
icolaou  editor,  SIAM  publications,  1986,  pp.  17-23. 

[12]  "Wave  propagation  and  heat  conduction  in  random  media",  in  Proceedings  of  International 
Congress  of  Mathematicians,  Berkeley  1986,  vol  2(1987)  pp.  1042-1051. 

[13]  "On  the  numerical  solution  of  a  nonlinear  stochastic  Helmholtz  equation  with  a  multigrid 
preprocessor"  (with  K.  Jordan  and  R.  Spigler)  Appl.  Math,  and  Computation  19  (1986)  pp.  145- 
157. 

[14]  "Focusing  and  Multifocusing  solutions  of  the  Nonlinear  Schroedinger  equation"  (with  B. 
LeMesurier,  C.  Sulem  and  P.  -L.  Sulem)  Physica  D  31  (May  1988)  78-102. 

[15]  "Pulse  reflection  by  random  media"  (with  Burridge,  Sheng  and  White)  in  "Non-Classical  Contin¬ 
uum  Mechanics"  J.  R.  Knopps  and  A.  A.  Laccy  editors,  Cambridge  University  Press,  1987,  pp.  3- 
21. 

[16]  "Waves  in  one-dimensional  random  media"  J.  L.  Lions  anniversary  volume,  in  Analyse  Mathcma- 
tique  et  Applications,  Gauthier-Villars,  Paris,  1988  pp  42 1-435. 
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[17]  "Probing  a  random  medium  with  a  pulse"  (with  Burridge,  Sheng  and  White)  SIAM  J.  Appl.  Math, 
49  (1989). 

[18]  "Wave  localization  characteristics  in  the  time  domain"  (with  Sheng,  White  and  Zhang)  Phys.  Rev. 
Letters  59  (Oct.  26  1987)  pp.  1918-1921. 

[19]  "Local  structure  of  the  self  focusing  singularity  of  the  nonlinear  Scroedinger  equation"  (with  B. 
LeMesurier,  P.  L.  Sulem  and  C.  Sulem)  Physica  D,  32  (1988)  pp.  210-226. 

[20]  "Rate  of  blow  up  for  solutions  of  the  nonlinear  Schroedinger  equation  at  critical  dimension"  (with 
M.  Landman,  P.  L.  Sulem  and  C.  Sulem)  Physical  Review  A,  38  (1988)  pp.  3837-3843. 

[21]  "Nonlinear  diffusion  limit  for  a  system  with  nearest  neighbor  interactions"  (with  M.  Z.  Guo  and  S. 
R.  S.  Varadhan)  Communications  of  Mathematical  Physics,  118  (1988)  pp.  31-59. 

[22]  "Frequency  content  of  randomly  scattered  signals  I”  (with  Asch,  Postel,  Sheng  and  White)  Wave 
Motion,  to  appear. 

[23]  "Frequency  content  of  randomly  scattered  signals  part  II:  Inversion"  (with  Postel,  Sheng  and 
White)  Wave  Motion,  to  appear. 

[24]  "Reflection  of  waves  generated  by  a  point  source  over  a  randomly  layered  medium"  (with  Kohler 
and  White)  Wave  Motion,  to  appear. 
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